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ABSTRACT 


A  multiscale  model  was  developed  to  link  the  hierarchies  of  human  bone  in  different 
length  scales.  Bone  has  a  unique  structure  displaying  large  stiffness  with  minimal  weight. 
This  is  achieved  through  a  hierarchy  of  complex  geometries  composed  of  only  three 
materials:  hydroxyapatite,  collagen  and  water. 

The  identifiable  structures  of  bone  are  hydroxyapatite,  tropocollagen,  fibrils, 
fibers,  lamellar  layers,  trabecular  bone,  cancellous  bone  and  cortical  bone.  A  spring 
model  was  used  to  evaluate  the  stiffness  of  collagen.  A  unit-cell  based  micromechanics 
model  analyzed  both  the  normal  and  shear  properties  of  fibrils,  fibers,  and  lamellar 
layers.  A  layered  composite  model  assessed  cortical  and  trabecular  bone  while  a  simple 
finite  element  model  was  used  to  evaluate  cancellous  bone. 

Modeling  bone  from  nanoscale  components  to  macroscale  structures  allows  the 
influence  of  each  structure  to  be  assessed.  It  was  found  that  the  distribution  of 
hydroxyapatite  within  the  tropocollagen  matrix  at  the  fibril  level  influences  the 
macroscale  properties  the  most.  Additionally,  the  model  allows  perturbations  to  the 
geometry  of  any  hierarchy  to  be  analyzed.  With  so  little  known  about  the  detailed 
structure  of  nanoscale  and  microscale  bone,  a  model  comprising  the  complete  hierarchy 
of  bone  can  be  used  to  help  validate  assumptions  or  hypotheses  about  structure. 
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I.  INTRODUCTION 


Biomaterials  are  complex  living  tissues  that  have  developed  through  evolutionary 
processes.  They  are  distinct  for  their  complex  hierarchies  built  by  simple  materials  [1]. 
There  are  a  limited  number  of  structural  materials  in  the  human  body,  but  living 
organisms  rely  on  composite  hierarchical  structures  to  achieve  macroscale  form  and 
function  [2],  Biomaterials  such  as  skin,  ligaments,  tendons,  muscles  and  bones  exhibit 
this  hierarchical  organization.  This  study  will  focus  on  the  structure  of  bone. 

Bone  is  the  main  structural  component  of  the  body.  It  supports  the  compressive 
loads  of  our  body  weight,  serves  as  lever  arms  to  transmit  muscle  contractions  to 
movement  and  provides  rigidity  for  each  limb  [3]— [6] .  Unlike  tendons  and  ligaments, 
bone  structures  support  compressive  loads,  as  well  as  bending,  torsional  and  shearing 
loads  [5],  Bone  structure  has  evolved,  adhering  to  the  biological  adage  that  form  follows 
function.  Bone’s  complex  hierarchies  create  a  load  bearing  structure  of  minimal  weight, 
capable  of  supporting  enormous  loads  [4].  Additionally,  bone  serves  as  the  storage  site 
for  bone  marrow,  calcium  and  phosphate  [3],  [5].  In  order  to  understand  the  structural 
properties  of  bone,  a  multi-scale  model  will  be  presented.  This  model  will  begin  at  the 
nano  level  and  continue  up  the  hierarchies  to  the  macroscale  level;  however,  for 
descriptive  purposes,  bone  will  be  introduced  from  the  macroscale  down. 

A.  MACROSCALE  BONE 

There  are  five  types  of  bones  in  the  human  body:  long,  short,  flat,  sesamoid  and 
irregular  [6].  Long,  short,  and  flat  bones  are  described  by  their  names.  Sesamoid  bones 
exist  within  a  tendon,  such  as  the  patella.  Irregular  bones,  examples  being  vertebrae,  fit 
into  none  of  the  other  categories  [6].  Bones  are  composed  of  two  distinct  types  of  bone: 
woven  bone  and  lamellar  bone  [3],  [7].  Woven  bone  is  a  disordered  array  of  collagen 
fibrils  [8].  Lamellar  bone,  which  is  organized  into  distinct  sheets,  is  an  ordered  array  of 
bone  fibers  in  a  matrix  material  [1],  It  is  the  dominant  form  of  bone  and  closely 
resembles  a  layered  fiber-reinforced  composite  [9]. 
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Lamellar  bone  can  be  further  divided  into  two  bone  structures:  cortical  bone  and 
cancellous  bone  [7],  [10].  Cortical  bone,  also  known  as  compact  bone,  is  a  dense  form  of 
bone.  It  is  composed  of  concentric  layers  of  lamellar  bone  with  varying  fiber  directions 
[3],  The  stacking  of  concentric  layers  forms  an  osteon,  or  Haversian  system,  and  is  the 
representative  subunit  of  cortical  bone  [5],  They  contain  a  Haversian  canal,  which 
supplies  the  blood  and  nerve  endings  to  bone,  at  their  center.  Additionally,  there  are  bone 
cells  named  osteoblasts  that  form  new  bone  and  matrix  material  [11],  These  are 
interconnected  by  a  network  of  canaliculi  that  help  transport  necessary  nutrients  and 
waste  materials  in  the  compact  bone  structures  [9].  These  osteons  are  densely  packed  in 
cortical  bone. 

Cancellous  bone,  otherwise  known  as  spongy  or  trabecular  bone,  is  a  very  porous 
material  [3],  [10].  Cancellous  bone  is  very  similar  in  form  to  metal  foam;  it  is  a  complex 
array  of  rods  and  thin  plates.  These  rods  are  known  as  trabeculae  and  are  composed  of 
lamellar  bone.  The  material  of  trabeculae  will  be  referred  to  as  trabecular  bone  while  the 
material  of  cancellous  bone  as  a  macroscale  material  will  be  referred  to  as  cancellous 
bone.  The  lamellar  bone  of  a  trabecula  has  a  uniform  fiber  direction  along  the  axis  of  the 
structure  [9].  This  is  caused  by  the  constant  absorption  and  re-growth  of  lamellar  material 
as  a  direct  result  of  the  stresses  applied  to  the  bone  [4],  [9], 

Cortical  bone  is  predominantly  found  at  the  cortex  of  bones  and  along  the  long 
axis  of  long  bones.  It  accounts  for  almost  80%  of  the  skeletal  system’s  weight. 
Cancellous  bone  is  present  along  the  inner  radius  of  long  bones,  within  the  heads  of  long 
bones,  and  acts  as  filler  under  the  cortex  of  cortical  bone.  Cancellous  bone  also  regulates 
the  metabolic  functions  executed  by  bone. 

B.  MICROSCALE  BONE 

The  lamellae  of  both  cortical  bone  and  trabecular  bone  are  composed  of  two 
microstructures:  fibrils  and  fibers  [3].  Bone  fibrils  are  an  extremely  ordered  array  of  bone 
mineral  and  tropocollagen  (TC)  molecules  [2],  [12]— [13].  The  bone  mineral  present 
within  fibrils  is  classified  as  intrafibrillar  mineral.  TC  molecules  form  a  lattice  structure 
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that  is  reinforced  by  the  bone  minerals  [13].  The  fibril  material  is  representative  of  a 
particle  reinforced  composite  [1]. 

Bone  fibers  are  an  ordered  array  of  mineralized  fibrils  [12],  The  individual  fibrils 
are  surrounded  by  a  coating  of  bone  mineral  known  as  extrafibrillar  mineral  [3],  [14], 
These  fibers  are  not  independent  structures,  but  exist  only  within  the  complex  lamellar 
bone  [13],  similar  to  a  fiber  reinforced  composite  [12],  Additionally,  the  fibers  are 
surrounded  by  an  extrafibrillar  matrix  that  contains  the  cells  and  non -structural  proteins 
necessary  for  metabolic  functions  [4], 

C.  NANOSCALE  BONE 

Bone  is  almost  entirely  composed  of  two  materials:  hydroxyapatite  (HA)  and 
collagen  [3],  HA,  also  known  as  bone  mineral,  is  a  small,  thin,  plate  structure  with  a  large 
elastic  modulus  [4],  [9].  The  collagen  species  of  bone  is  95%  collagen  I  and  a  5% 
assortment  of  other  collagen  species.  The  main  constituent  of  collagen  I  are  TC 
molecules.  TC  molecules  are  a  triple  helix  molecule  [8]  composed  of  three  polyproline  II 
helices  [1 1]. 

D.  HIERARCHICAL  OVERVIEW 

From  just  three  basic  materials,  HA,  TC,  and  water,  a  complex  hierarchy  of  bone 
is  built.  This  structure  has  evolved,  and  continues  to  evolve,  to  support  the  stresses  it 
encounters  [4],  The  astounding  aspect  of  bone  is  that  its  material  properties  have  been 
produced  through  the  complex  geometry  of  simple  materials  [9].  This  hierarchy  is  broken 
down  and  visualized  in  Figure  1 . 
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Nano  scale 


(a)  Tropocollagen 


(b)  Hydroxyapatite 


Microscale 


Macroscale 


(f)  Trabecular  Bone 


(g)  Cortical  Bone 


Figure  1 .  Bone  hierarchies  as  classified  into  nanoscale  (a-b),  microscale  (c-d), 
and  macroscale  (e-g)  structures,  (a)  Tropocollagen  molecules  comprised 
of  three  polyproline  II  helices,  (b)  Hydroxyapatite  crystals,  (c) 
Organized  fibril  lattice  structures,  (d)  Mineralized  fibers,  (e)  Lamellar 
bone  constructed  of  layered  fibrils  in  a  matrix,  (f)  Trabecular  bone  with 
parallel  fiber  layers,  (g)  Cortical  bone  with  alternating  fiber  orientations 
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II.  NANOSCALE 


The  main  structural  nano  components  of  bone  are  hydroxyapatite  and 
tropocollagen.  These  two  materials  serve  as  the  building  blocks  for  the  subsequent  layers. 
In  addition  to  hydroxyapatite  (HA)  and  tropocollagen  (TC)  there  are  microscopic 
inclusions  of  non-collagenous  proteins  (NCPs)  and  proteoglycans.  The  NCPs  serve 
metabolic  functions  and,  along  with  the  proteoglycans,  act  as  a  matrix  fdler  [9],  Since  the 
NCPs  and  proteoglycans  serve  no  mechanical  functions  and  comprise  a  small  volume 
percentage,  they  are  ignored  in  this  study.  The  following  chapter  will  describe  in  detail 
HA  and  TC,  and  will  propose  a  mechanical  model  for  the  TC  structure. 

A.  HYDROXYAPATITE 

HA  is  a  crystalline  solid  composed  of  calcium  phosphate.  It  is  the  only  ceramic 
material  created  naturally  within  the  body  [8],  [15].  HA  is  sometimes  referred  to  as  the 
mineral  phase,  or  bone  mineral,  and  takes  a  polycrystalline  form  as  small  thin  plates 
within  the  body  [4],  [16],  HA  is  the  direct  cause  of  the  stiffness  associated  with  bone. 

1.  Chemistry 

The  chemical  formula  of  HA  is  Caio(P04)6(OH)2  [8],  The  growth  of  HA  is 
regulated  by  heterogeneous  nucleation  factors,  which  are  thought  to  be  concentrated 
within  the  gap  zones  [4],  These  nucleation  factors  promote  the  growth  of  the  mineral 
phase  and  contribute  to  the  highly  ordered  structure  of  the  micro-level  hierarchy  of  bone. 
Additionally,  HA  exhibits  a  hexagonal  close-packed  (HCP)  crystal  lattice.  A  single 
crystal  can  be  grown  in  a  lab  and  tested  to  determine  material  properties  [16].  The  HCP 
crystal  lattice  of  HA  is  shown  in  Figure  2. 
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Figure  2.  Hydroxyapatite  crystal  lattice  structure:  hexagonal  close-packed. 

2.  Material  Properties 

From  nanoindentation  testing  completed  on  single  HA  crystals,  the  stiffness  in  the 
[0001]  direction  was  found  to  be  greater  than  in  the  [1010]  direction  [16].  The  results  of 
these  tests  are  shown  in  Table  1 ,  along  with  several  other  experimental  results. 


Table  1 .  Hydroxyapatite  material  values. 


Direction 

E  [GPa] 

V 

Source 

[0001] 

150.38 

Nanoindentation  [16] 

[1010] 

143.56 

Nanoindentation  [16] 

114 

0.27 

Ultrasonic  [17] 

165 

[8] 

The  results  from  [17]  are  representative  of  ultrasonic  moduli  testing  completed  on 
pine  powdered  materials.  These  results  assumed  an  isotropic  material  in  pure 
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polycrystalline  form.  The  estimates  from  [8]  cite  no  methods  or  sources.  The  values  from 
the  nanoindentation  experiments  will  be  used  in  this  study. 

Furthermore,  the  HCP  structures  are  organized  into  small  plates  ranging  from 
1.5-5  nni  thick  [9].  Although  the  size  of  FI  A  crystals  vary  based  on  location,  mineral 
density  and  time  allowed  for  growth,  the  average  crystal  size  is  widely  agreed  upon  at  50 
nm  x  25  nm  x  3  nm  [1],  [3],  [4],  [8],  [16],  [18],  [19]. 

3.  Non-structural  Function 

In  addition  to  providing  structural  support  to  the  skeletal  system,  HA  serves  as  the 
storage  site  for  ions  [16].  HA  is  responsible  for  99%,  90%,  90%  and  50%  of  the  body’s 
store  of  calcium,  phosphorus,  sodium  and  magnesium,  respectively  [4],  HA  has  a  large 
surface  area  to  volume  ratio,  allowing  for  rapid  absorption  and  dissolution  when  ions  are 
needed  [4],  [16].  This  relationship  is  important  in  determining  effects  of  bone  diseases 
and  disorders  that  affect  bone  density  and  calcium  absorption. 

B.  TROPOCOLLAGEN 

The  tropocollagen  molecule  is  the  basic  molecular  unit  of  collagen.  There  are  28 
different  variations  of  collagen;  each  is  used  for  a  different  physiological  purpose  within 
the  body  [20].  Collagen  I  is  the  main  collagen  species  found  in  bone,  comprising  95%  of 
the  total  collagen  [21], 

1.  Tropocollagen  Structure 

Collagen  I  is  composed  of  three  helical  protein  strands:  two  alpha-1  type-1  strands 
and  one  alpha-1  type-2  strand  [20],  [22]-[26],  Each  strand  differs  in  its  amino  acid 
sequence,  but  each  maintains  a  similar  helical  structure.  Both  alpha-1  type-1  and  alpha-1 
type-2  are  left  handed  polyproline  II  helices,  are  300  nm  long,  and  contain  approximately 
1000  amino  acids  [11],  [22],  [26]-[32],  Additionally,  the  polyproline-II  (PPII)  helices  of 
the  alpha- 1  type-1  and  alpha- 1  type-2  collagen  contain  a  common  repeating  subunit:  Gly- 
X-Y  [20],  [23]-[25],  [27],  [28],  [33].  Gly  is  the  amino  acid  glycine,  the  smallest  of  the 
amino  acids.  With  this  repeat  for  each  strand,  a  glycine  amino  acid  is  along  the  central 
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axis  of  the  tropocollagen  molecule  [23],  [34],  [35].  This  allows  the  three  strands  to  orient 
themselves  in  such  a  way  as  to  form  stabilizing  hydrogen  bonds. 

The  combination  of  the  alpha- 1  type-1  and  alpha- 1  type-2  helices  creates  the  right 
handed  helix  of  the  TC  molecule.  Additionally,  a  pattern  exists  between  the  X  and  Y 
position  of  the  three  amino  acid  subunit.  Two  common  repetitions  are  Gly-Pro-Y  and 
Gly-X-Hyp,  where  Pro  and  Hyp  represent  the  amino  acids  proline  and  hydroxyproline, 
respectively  [23]-[25],  [30],  A  common  model  for  the  TC  molecule  utilizes  the  repeating 
sequence  Gly-Pro-Hyp,  which  is  shown  in  Figure  3. 


Figure  3.  Triple  helix  of  collagen  with  a  Gly-Pro-Hyp  model;  colored  ribbons 

indicate  individual  PPII  helix,  from  [25]. 
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2.  Spring  Model 

In  order  to  discern  the  stiffness  of  a  single  molecule  of  tropocollagen,  a  basic 
mechanic  model  of  a  helical  spring  was  used.  This  model  can  be  used  because  each  left- 
handed  polyproline  helix  of  the  TC  molecule  can  be  independently  treated  as  a  spring. 
The  stabilizing  hydrogen  bonds  and  van  der  Waals  attraction  between  the  three  helices 
helps  to  prevent  buckling  and  to  ensure  that  all  three  helices  deform  symmetrically. 
Additional  protection  against  buckling  is  provided  by  the  close  packing  of  the  TC 
molecules  within  collagen  fibrils,  to  be  discussed  later.  The  stiffness  of  the  TC  molecule 
is  the  summation  of  the  stiffness  of  the  alpha- 1  type-1  and  alpha- 1  type-2  strands. 
Equation  1  was  used  to  find  the  spring  stiffness  of  each  helix. 


F  _  Gd 4 
y  "  8 D3NA 


(1) 


where  k  represents  the  spring  rate,  G  represents  the  shear  modulus,  d  represents  the  wire 
diameter,  D  represents  the  mean  coil  diameter,  and  NA  represents  the  number  of  coils. 


3.  Model  Derivations 

To  determine  the  quantitative  values  for  each  protein  helix,  it  was  assumed  that 
each  repeating  subunit  is  represented  by  Gly-Pro-Hyp.  This  representative  subunit 
provides  a  relatively  stable  position  for  each  amino  acid  and  allows  for  proper  cross- 
linking  and  hydrogen  bonds  to  form  between  strands  [25],  [27],  [30],  [36].  Additionally, 
it  is  the  presence  of  the  amino  acids  proline  and  hydroxyproline  that  most  greatly 
influence  the  stiffness  of  collagen  [37],  Although  in  reality  the  combination  of  proline 
and  hydroxyproline  comprise  only  one-sixth  of  the  amino  acid  content,  the  chemical  and 
molecular  level  effects  of  the  excess  proline  and  hydroxyproline  can  be  ignored  since  this 
model  relies  predominantly  on  the  structure  and  not  the  molecular  interactions  [34], 
Additionally,  modeling  each  subunit  with  this  sequence  matches  the  correct  bond  angles 
and  axial  repeat,  10/3,  for  the  majority  of  the  polyproline  helices  [25].  Recent  studies 
have  debated  the  predominance  of  an  axial  repeat  of  7/2  over  portions  of  the  PPII; 
however,  the  implications  of  this  variation  are  more  important  for  molecular  dynamics 
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studies  where  strain  energies  and  steric  effects  of  the  aniino  acid  interactions  are 
computed  and  analyzed  [34],  [35], 

The  wire  diameter  was  taken  to  be  the  average  spacing  between  residues,  0.286 
mil  [26],  [34],  [35],  The  coil  diameter  was  derived  from  the  diameter  of  a  TC  molecule.  A 
TC  molecule  is  1.5  mn  in  diameter,  thus,  a  single  PPII  helix  was  assumed  to  have  a  mean 
coil  diameter  of  0.5  run  [25],  [27],  [28],  The  number  of  coils  was  derived  from  a  constant 
axial  repeat  of  10/3  [23],  [30],  With  1000  amino  acids,  this  yields  300  coils. 

The  shear  modulus  was  the  most  difficult  assumption  to  ascertain,  since  no  data  is 
available  as  to  the  shear  modulus  of  a  single  chain  amino  acid  helix.  In  order  to  create  a 
valid  estimate,  the  bond  energy  of  the  backbone  of  a  single  subunit  was  calculated. 
Within  the  Gly-Pro-Hyp  triplet,  there  are  six  C-N  bonds  and  three  C-C  bonds  along  the 

backbone.  With  bond  energies  of  5.06  x  10~19  ^  and  5.75  x  10~19 — — — ,  respectively, 

bond  bond 

the  total  energy  of  the  backbone  is  then  4.76  x  10-18 — — — .  With  a  subunit  volume  of 

subunit 

7.02  x  10-29  m3,  the  shear  stress  of  the  backbone  is  computed  as  energy  over  volume, 
which  equates  to  67.9  GPa.  This  is  close  to  the  typical  values  of  aluminum  (70  GPa),  and 
thus  not  unreasonable.  The  calculation  above  neglects  all  strengthening  effects  of  cross- 
linking.  non-backbone  atoms  and  other  atomistic  considerations,  but  still  provides  a  basic 
starting  point  to  calculate  the  overall  stiffness  of  collagen. 

The  resulting  values  computed  for  Equation  1  are  displayed  in  Table  2. 

Table  2.  Variable  values  for  Equation  1 . 


Variable 

Value 

Units 

G 

67.9 

GPa 

d 

0.286 

mn 

D 

0.5 

run 

300 

k 

0.0015 

N 

m 

10 


4. 


Model  Results 


The  stiffness  of  each  helix  can  then  be  found  by  equating  Young’s  Modulus  to 
spring  rate. 


E 


(2) 


From  the  values  stated  in  Table  2,  a  spring  stiffness  of  0.0015  — was  calculated 

m 

from  Equation  1.  The  length  and  diameter  of  a  PPII  helix  are  300  nm  and  0.5  nm, 
respectively.  However,  due  to  the  helical  twist  shown  in  Figure  3,  the  cross-sectional 
diameter  is  defined  as  0.7  nm.  When  the  length  of  300  nm  and  cross  sectional  area  of 
3.85  x  10~19  nr  are  applied  to  Equation  2,  a  stiffness  of  1.18  GPa  is  found  for  each 
protein  helix.  Since  each  TC  molecule  is  the  sum  of  three  helices,  the  stiffness  of  a  TC 
molecule  is  3.54  GPa. 

A  comparison  of  the  results  of  this  study  with  known  values  from  both 
experimental  data  and  theoretical  analysis  are  shown  in  Table  3. 
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Table  3.  Comparison  of  values  for  tropocollagen  stiffness. 


Stiffness 

[GPa] 

Method 

Source 

3.54 

Spring  model 

[This  Study] 

1.2 

Property  used  in  finite  element  model 

[1] 

2.8 

Property  used  in  finite  element  model 

[11] 

6-16  (Average  6) 

Molecular  dynamics 

[38] 

7 

Atomistic  modeling 

[39] 

2.4 

Atomistic  modeling 

[40] 

3 

X-ray  diffraction 

[41] 

9 

Brillouin  light  scattering 

[42] 

5.1 

Brillouirr  fight  scattering 

[43] 

3 

Estimate  from  persistence  length 

[44] 

0.35-12 

Estimate  from  persistence  length 

[40],  [45] 

4.8±1 

Molecular  dynamics 

[46] 

1.4 

Debye-Waller  factor 

[30] 

(PHG)1  11.37 

Molecular  dynamics 

[36] 

(PPG)2  13.43 

Molecular  dynamics 

[36] 

'PHG  is  a  Gly-Pro-Hyp  model 
2PPG  is  a  Gly-Pro-Pro  model 


These  results  demonstrate  that  the  helical  spring  model  produces  a  valid  answer 
for  TC  stiffness,  thus  confirming  the  accuracy  of  this  method.  What  this  model  does  not 
provide  is  the  transverse  properties  or  a  Poisson’s  ratio  for  TC.  Thus,  we  must  assume 
that  collagen  at  upper  hierarchies  is  isotropic  with  a  Poisson’s  ratio  of  0.35  [1],  This 
relatively  high  Poisson’s  ratio  is  due  to  the  presence  of  water,  which  serves  as  a 
stabilizing  agent  along  the  length  of  the  helix  [1 1],  [33].  [37]. 

C.  NANOSCALE  CONSIDERATIONS 

The  nanoscale  model  presented  is  a  mechanical  approach  to  a  traditionally 
chemical  or  biological  problem,  as  modeling  of  molecules  is  generally  restricted  to 
molecular  dynamics  or  atomistic  calculations.  However,  with  an  ordered  hehcal  structure 
constrained  on  all  sides  to  prevent  buckling,  a  simple  mechanical  model  produces 
accurate  results.  This  model  neglects  the  mechanical  properties  of  water  bound  within  the 
helices. 
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III.  MICROSCALE 


The  previous  chapter  defined  the  material  properties  of  the  elemental  composite 
materials  of  bone.  From  those  two  simple  constituents,  a  more  complex  hierarchy  is 
produced.  The  microscale  structures  of  the  ensuing  hierarchy  are  highly  organized  and 
highly  regular  structures  of  collagen  and  hydroxyapatite.  In  order  to  calculate  the 
representative  material  properties  of  these  microstructures  a  micromechanics  model  is 
utilized. 

A.  BONE  FIBRIL 

The  fibril  is  the  smallest  of  the  microstructures.  It  is  represented  by  a  staggered 
array  of  hydroxyapatite  and  tropocollagen  molecules.  These  staggered  arrays  are 
dependent  on  the  binding  and  cross-linking  tendencies  of  tropocollagen  molecules. 

1.  Two-Dimensional  Fibril  Structure 

The  regular  organization  of  fibrils  is  caused  and  propagated  during  the  expansion 
and  re-growth  of  fibrils,  known  as  fibrillogenesis.  When  a  collagen  molecule  is  produced 
it  has  a  C-terminal  and  N-terminal,  so  named  for  the  peptides  cleaved  during  the 
transition  from  procollagen  to  collagen  [26],  [32],  The  amino  acid  sequence  at  each 
terminal  of  the  TC  defines  a  polar  reference  for  each  molecule,  as  each  end  exhibits 
different  cross-linking  tendencies.  When  the  collagen  molecules  are  aligned  near  each 
other  the  C  and  N  ends  overlap  to  produce  an  organized  structure.  This  preferential 
alignment  is  not  unique  to  the  ends  of  the  TC  molecules;  the  entire  length  of  the  TC 
molecule  contains  segments  of  amino  acids  preferential  to  cross-linking  with  adjacent 
tropocollagen  molecules  [33].  These  preferential  segments  form  a  collagen  network  with 
a  40  mil  gap  and  a  regular  periodicity  of  67  nm  [29],  [32],  [33],  [38],  [47]-[50],  These 
gap  regions  are  preferential  to  the  growth  of  HA,  allowing  for  the  inclusion  of  the 
reinforcing  particles  within  the  fibril  composite.  This  structure  is  shown  in  Figure  4. 
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Figure  4.  Visualization  of  67  nm  periodicity  in  bone  fibril.  Dark  regions 
denote  hydroxyapatite  crystals,  light  regions  depict  tropocollagen. 

Figure  4  is  not  to  scale.  The  height  of  the  crystals  and  TC  matrix  are  3  nm,  while 
the  F1A  crystals  are  50  mu  in  length.  This  image  is  used  to  portray  the  periodic  packing 
that  is  known  to  occur.  This  framework  provides  a  regular  structure  driven  by  the 
molecular  sequence  of  the  tropocollagen  molecules.  Additionally,  this  cross-linking 
provides  a  stabilizing  effect  on  the  tropocollagen  phase  of  the  fibril  [29],  [37].  The 
presence  of  cross-linking  prevents  the  individual  molecules  from  shearing.  This  in  turn 
increases  strength  and  brittleness  of  the  structure  [11],  [29].  Furthermore,  the  presence  of 
cross  linking  allows  the  collagen  phase  to  be  represented  as  a  single,  continuous,  material 
in  a  simplified  model. 

An  important  quality  of  the  tropocollagen  and  hydroxyapatite  microstractures  is 
that  the  HA  crystals  are  aligned  parallel  to  the  microfibril  direction.  This  is  due  to  the 
growth  behavior  of  the  mineral  phase  within  the  gap  zones. 

a.  Heterogeneous  Nucleation  Factors 

Heterogeneous  nucleation  factors  assist  in  the  initial  deposition  of  HA  at  a 
molecular  level  [4].  An  important  nucleation  factor  is  bone  sialoprotein  (BSP).  BSP  binds 
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directly  to  the  ends  of  collagen  I  molecules  and  can  be  found  in  high  concentrations 
within  the  gap  zones  of  microfibrils  [51].  BSP  has  been  found  to  be  a  nucleator  of 
hydroxyapatite  formation,  creating  a  HA  crystal  within  a  supersaturated  medium  [51]. 

b.  Hydroxyapatite  Growth 

Hydroxyapatite  crystals  have  shown  preferential  growth  orientation  related  to  the 
substrate  chemistry  on  which  they  are  grown  [52],  Sato  et  al.  [52]  believe  that  the  growth 
orientation  was  related  to  the  alignment  of  carboxyl  groups  within  an  organic  film.  Chen 
et  al.  [53]  further  showed  that  the  hydrothermal  crystallization  of  hydroxyapatite  can 
produce  an  oriented  attachment  through  the  use  of  ordered  organic  materials. 

This  process  is  similar  to  the  growth  of  crystals  in  vivo.  The  gap  zones  are  first 
supersaturated  with  mineral  components.  Then,  the  BSP  proteins  initiate  crystal  growth 
on  the  ends  of  the  collagen  I  molecules.  Subsequent  crystal  growth  occurs  in  an  ordered 
manner  that  aligns  all  crystals  in  parallel  with  the  microfibril  axis  [51].  Since  all  crystals 
are  grown  in  the  same  orientation,  the  stiffness  values  shown  in  Table  1  can  be  used  as 
material  properties  for  a  transversely  isotropic  material. 

One  issue  arises  with  the  physical  models  presented  based  on  TC  cross-linking: 
the  40  nm  gaps  present  in  the  collagen  network  are  smaller  than  the  50  nm  x  25  nm  x 
3  nm  HA  crystals.  Researchers  believe  that  the  HA  must  surround  the  ends  of  the  TC 
molecules  and  grow  between  adjacent  molecules  [50],  [54],  This  helps  to  prevent 
shearing  of  the  HA-TC  bonded  surfaces.  The  width  of  the  crystals  represents  the  next 
challenge  in  modeling  and  determining  the  3D  structure  of  collagen  fibrils. 

2.  Three-Dimensional  Fibril  Structure 

The  regular  staggered  array  of  a  collagen  fibril  is  easily  visualized  in  two 
dimensions,  however,  three  dimensional  fibrillogenesis  in  vivo  is  still  a  debated  topic 
[32],  [33].  Electron  micrograph  images  of  tendon  fibrils  have  shown  different  three- 
dimensional  structures,  as  seen  in  Figure  5.  Tendon  fibrils  are  very  similar  in  structure  to 
bone  fibrils,  allowing  for  a  feasible  interpretation  of  HA  packing  models. 
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Figure  5.  Electron  micrograph  images  of  3D  structures  of  collagen  fibril,  from 
[26].  (a)  Linear  crystal  pattern,  (b)  Twisting  crystal  pattern. 

For  many  years  a  linear  fibril  orientation  has  been  the  accepted  model  [19],  [32], 
[33],  [55],  [56].  This  model,  shown  in  Figure  5a,  was  viewed  as  a  long  thin  filament  with 
alternating  bands  of  mineral  rich  phase,  sometimes  referred  to  as  the  gap,  and  mineral 
deficient  phase,  sometimes  referred  to  as  the  overlap.  The  mineral  rich  phases  have 
approximately  twice  the  mineral  content  of  the  mineral  deficient  phase  [55].  Figure  6 
helps  visualize  the  three-dimensional  packing  of  HA  crystals  in  the  linear  fibril  model. 
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Although  the  linear  model  presents  a  viable  solution  to  the  proposed  theory  of 
three-dimensional  fibrillogenesis,  there  has  not  been  any  evidence  proving  lateral  growth 
is  a  purely  linear  process.  A  recent  study  conducted  on  the  in  vitro  organization  of  fibril 
growth  discovered  that  fibrils  grow  laterally  in  4  nm  steps  [47].  Each  lateral  step  is 
defined  by  electrostatic  attraction  during  the  formation  of  the  microscale  fibril  [47].  The 
study,  however,  did  not  define  lateral  growth  as  it  relates  to  longitudinal  growth.  A  more 
recent  model  being  considered  is  that  of  a  twisting  crystalline  structure  [9].  A  spiral 
provides  for  a  more  even  distribution  of  HA  crystals.  This  feature  is  shown  in  Figure  5b. 
In  the  model  presented,  it  is  assumed  that  the  67  nm  periodicity  that  determined  the  two- 
dimensional  stacking  also  directs  the  three-dimensional  pattern.  For  this  reason,  the 
lateral  periodicity  will  also  be  67nm.  A  simplified  model  is  shown  in  Figure  7. 
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Figure  7.  Twisting  fibrillar  model. 


The  two  proposed  fibril  models  shown  in  Figure  6  and  Figure  7  are  meant  to  help 
portray  a  simplified  visualization  of  F1A  crystal  distribution.  The  organization  of  bone 
fibril  in  vivo  is  a  much  more  complex  occurrence  that  includes  non-collagenous  proteins, 
errors  in  collagen  stacking,  errors  in  amino  acid  sequencing  and  continuous  absorption 
and  re-growth  of  fibrils.  Moreover,  biological  materials  are  never  pure  substances.  The 
collagen  used  in  bone  is  only  95%  collagen  I  [21].  Piez  and  Miller  [26]  postulate  that  the 
other  collagen  species  help  to  regulate  the  fibrillogenesis  process,  but  they  may  also 
introduce  irregularities  into  the  structure. 

Fibrils  are  generally  accepted  to  be  circular,  with  a  diameter  ranging  from  50  nm 
to  300  nm  [57],  [58].  It  is  possible,  in  these  small  diameter  fibrils,  that  only  one  crystal 
width  is  present  across  the  structure.  This  would  validate  the  use  of  the  linear  fibril 
model.  However,  in  larger  fibrils  the  diameter  allows  for  either  fibril  structure  to  be 
present.  In  order  to  explore  the  two  models  discussed,  a  micromechanical  model  will  be 
modified  to  characterize  both  assemblies.  The  results  will  then  be  compared  to 
experimental  and  theoretical  data  to  determine  validity. 
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3. 


Micromechanical  Model 


Kwon  and  Kim  [59]  developed  a  micromechanics  model  for  the  analysis  of 
composite  structures.  This  model  is  comprised  of  a  rectangular  prism  divided  into  8 
quadrants,  or  subcells,  shown  in  Figure  8. 


Figure  8.  Unit  cell  of  micromechanics  model. 

This  unit  cell  represents  a  composite  material,  where  each  subcell  can  be 
identified  as  a  different  material.  This  model  can  solve  for  equivalent  composite  stiffness, 
rigidity,  Poisson’s  ratio  and  coefficient  of  thermal  expansion. 

a.  Model  Derivation 

In  Figure  8,  the  numbers  on  the  unit  cell  faces  indicate  the  subcell  numbers  while 
a—c  indicate  subcell  dimensions.  Additionally,  each  subcell  is  assumed  to  have  uniform 
stress  and  strain  so  that  the  equilibrium  at  each  interface  can  be  defined  [59]— [6 1  ] .  The 
following  equations  define  the  equilibrium  stress  state  at  the  interfaces  of  the  subcells. 
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The  subscripts  denote  stress  components  according  to  the  axis  in  Figure  8  while 
the  superscripts  indicate  subcell  number.  In  addition  to  stress  equilibrium,  deformation 
compatibility  must  be  met  for  the  unit  cell.  Deformation  compatibility  defines  the 
following  set  of  equations. 
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Furthermore,  each  subcell  has  a  constitutive  equation  defining  strain  [61]. 


n  n  n  _ _ n  \ 

ij  ijkl  kl  ij 


(7) 


The  superscripts  represent  the  subcell  number  and  subscripts  represent  tensor  notation. 
C*Jkl  is  the  stiffness  tensor  and  a,"  is  the  thermal  expansion  tensor.  For  this  model,  the 

thermal  expansion  of  the  materials  is  ignored.  This  can  be  assumed  because  the  internal 
temperature  of  the  human  body  is  highly  regulated.  This  internal  temperature  regulation 
results  in  a  negligible  thermal  stress  component. 


20 


The  volume  fraction  of  each  subcell  can  be  calculated  from  the  dimensions  of  the 
unit  cell  and  represented  in  vector  form  as  V"  or  matrix  form  as  [V].  The  total  unit  cell 

stress  and  strain  can  then  be  found  by  averaging  the  subcell  stresses  and  strains  based  on 
subcell  volume  [59]— [6 1  ]. 
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The  overbars  in  Equations  8  and  9  represent  the  unit  cell  values.  With  Equations 
3,  5,  and  9,  a  matrix  equation  can  be  formed  for  the  normal  properties  of  the  composite, 
as  shown  in  Equation  10. 

bM={/}  do) 

[T]  is  a  24x24  matrix  with  the  combination  of  three  vectors  containing  stress 
relationships,  strain  relationships  and  constitutive  equations. 

[rf=[[n]  [T2]  [re]] 

[ri]  is  a  24x12  matrix  of  normal  stress  relationships,  \T 2]  is  a  24x9  matrix  of  normal 
strain  relationships  and  [723]  is  a  24x3  matrix  of  constitutive  equations.  {/}  is  a  24x1 

column  vector  composed  of  a  21  x  1  column  containing  zeros,  representing  the  stress  and 
strain  equilibrium,  and  a  3  x  1  column  containing  the  effective  normal  strains  [59]. 

{/}  *  =  [0]  ^£u  £ 22  £33 

The  assumption  that  thermal  stresses  are  zero  has  greatly  simplified  the  model 
from  [59].  Algebraic  arrangement  of  [r]  within  Equation  12  produces  the  following 
relationship. 

m  (id 
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The  inverse  of  [r]can  be  further  broken  down  into  three  matrices.  These  can  then 
be  substituted  into  Equation  1 1  to  obtain  the  following. 

PIMM  [“1  [®D 

_!  (12) 

{^}=[i?3]^ll  S  22  £33] 

Furthermore,  Equation  8  can  be  written  as: 

<13> 

With  Equations  7  and  12  substituted  into  Equation  13,  the  following  unit  cell 
stress  equation  is  found. 

<14) 

[£]  is  a  matrix  of  the  inverse  of  the  subcell  compliance  tensors  stated  in  Equation  7. 
From  Equation  14  the  unit  cell  stiffness  can  be  calculated. 

£]  =  [F][£][7?3]  (15) 

The  E  matrix  found  through  Equation  15  is  the  3x3  matrix  of  the  unit  cell  stiffness. 

With  these  values  the  Poisson’s  ratio  of  the  unit  cell  can  also  be  found.  A  similar  method 
can  be  utilized  to  solve  for  the  shear  properties  of  the  composite  material.  However, 
different  stress  relationships  and  deformation  values  are  utilized.  Equation  1 0  remains  the 
base  equation,  but  with  a  slight  variation.  The  shear  properties  will  be  defined  through 
Equation  16. 

[TT]{s}  =  {ff}  (16) 

The  [TT ]  matrix  utilized  for  shear  properties  is  composed  of  three  sub-matrices, 
as  with  normal  deformation. 

[7Y]r  =  [[m]  [7T2]  [2T3]] 
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[JTl]  is  a  24x18  matrix  of  shear  stress  relationships,  [JT2]  is  a  24x3  matrix  of  shear 

strain  relationships  and  [7T3]is  a  24x3  matrix  of  shear  constitutive  equations.  { }  is  a 

24x1  column  vector  composed  of  a  21xl  column  containing  zeros,  representing  the 
stress  and  strain  equilibrium,  and  a  3x1  column  containing  the  effective  shear  strains 
[59]. 


iff) 


-1 


[0]  [/n  r  13  y 23 


Derivation  of  composite  properties  follows  the  process  shown  through  Equations 
11-14,  however,  shear  calculations  utilize  the  [7Y]  and  {ff}  matrices.  The  resulting 
equation  solves  for  the  composite  shear  properties. 

'G]  =  [E][G][^3]  (17) 


is  representative  of  the  inverse  of  the  compliance  matrix  for  the  composite.  [iTK3]  is 


derived  from  the  inverse  of  the  [ TT ]  matrix.  This  model  assumes  independent  normal 

and  shear  properties.  Modification  is  simple  because  the  only  inputs  required  are  the 
dimensions  of  the  unit  cell  and  the  compliance  tensor  of  each  material  used.  The  coding 
of  the  model  allows  for  each  cell  to  have  a  different  material  property,  if  needed. 


b.  Subunit  Quantification 

For  the  purpose  of  this  work,  the  collagen  network  is  assumed  to  stagger  its  array 
with  at  least  a  25  nm  width,  so  as  to  allow  for  the  generally  accepted  dimensions  of  the 
HA  crystals.  However,  the  micromechanical  model  will  provide  relevant  results  for  both 
the  lateral  packing  model  and  the  twisted  packing  model  mentioned  above.  For  both 
models,  the  material  properties  shown  in  Table  4  will  be  used. 

In  addition  to  using  two  different  models  to  explore  the  crystal  arrays  within 
fibrils,  the  presence  of  water  will  be  taken  into  account  with  the  calculation  of  material 
properties.  The  presence  of  water  produces  a  bi-modulus  composite  material,  as  water  is 
only  considered  in  compression.  While  the  volume  percent  water  is  much  higher  in 
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tendons  and  ligaments,  bone  is  known  to  be  approximately  10-25%  water  [11],  Some  of 
this  is  thought  to  be  found  within  the  nanoscale  tropocollagen  and  hydroxyapatite.  Water 
serves  as  a  binding  and  stabilizing  agent  within  the  triple  helix  of  tropocollagen  [33], 
[37].  Additionally,  small  amounts  of  water  are  tightly  bound  within  the  HA  crystal  [11], 
However,  the  remaining  water  is  assumed  be  held  within  the  various  hierarchies.  For  the 
pmpose  of  the  fibrillar  model,  each  unit  cell  is  assumed  to  have  8%  volume  of  water. 
This  is  calculated  by  adjusting  the  stiffness  of  collagen  to  account  for  water  in 
compression  and  void  space  in  tension.  Tropocollagen  in  compression  and  tension 
exhibits  the  different  properties  shown  in  Table  4. 


Table  4.  Fibril  model  component  material  properties. 


Material 

E!  [GPal 

E23  [GPal 

V 

G  [GPal 

Hydroxyapatite 

150.38 

143.56 

0.23 

59.744 

Tropocollagen  (compression) 

3.428 

3.428 

0.35 

1.270 

Tropocollagen  (tension) 

3.256 

3.256 

0.35 

1.206 

(1)  Linear  Fibril  Subunit 

The  linear  packing  model  is  a  lateral  repeat  of  the  two-dimensional  fibril  array 
shown  in  Figure  4.  Each  ciystal  shown  in  the  model  has  a  length  of  50  mn,  a  width  of 
25  mn.  and  a  height  of  3  mn.  The  cross  section  used  to  create  a  repeating  subunit  of  the 
array  is  shown  in  Figure  9. 
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Figure  9.  Subunit  applied  to  the  linear  fibril  model. 


This  small  section  is  repeated  throughout  the  array  in  67  nm  increments.  The 
dimensions  as  they  apply  to  the  unit  cell  are  listed  in  Table  5. 


Table  5.  Unit  cell  dimension  for  linear  fibril  model. 


Dimension 

Value  [mn| 

ai 

50 

°2 

17 

\ 

25 

b2 

3 

ci 

3 

C2 

9 

This  model  allows  for  the  characteristic  67  nm  periodicity  in  the  1  direction,  a 
3  nm  buffer  between  each  subsequent  HA  crystal  in  the  2  direction,  and  the  vertical 
spacing  in  the  3  direction  associated  with  the  linear  packing  model.  Subcell  1  is  assigned 
HA  properties;  the  remaining  subcells  are  assigned  TC  properties. 


(2)  Twisting  Fibril  Subunit 
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The  twisting  fibril  takes  into  account  a  67  urn  periodicity  in  the  2  direction  as  well 
as  the  3  direction.  Additionally,  because  each  repeated  stack  in  the  2  direction  contains 
the  same  67  tun  periodicity  in  the  3  direction,  the  subunit  is  shortened  in  the  3  direction. 
Figure  10  shows  the  representative  area  used  for  the  twisted  fibril  subunit. 


Figure  10.  Subunit  employed  in  the  twisting  fibril  model. 


The  dimensions  associated  with  the  twisting  model  are  shown  in  Table  6.  It 
should  be  noted  that  the  dimensions  for  the  nl ,  a2,  bx ,  and  c,  lengths  are  restrained  due 
to  the  67  run  periodicity  in  the  1  dir  ection  and  the  known  size  of  the  HA  crystals. 

Table  6.  Unit  cell  dimension  for  twisting  fibril  model. 


Dimension 

Value  [mn] 

ai 

50 

°2 

17 

25 

K 

75 

e? 

3 

C2 

6 
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For  the  twisting  model,  subcell  1  is  assigned  HA  properties  and  the  remaining 
subcells  are  assigned  TC  properties. 

c.  Model  Results 

The  micromechanics  model  calculates  the  orthotropic  properties  that  result  from  a 
non-symmetric  unit  cell.  However,  in  vivo  we  can  expect  to  find  fibrils  with  transverse 
isotropic  properties.  This  is  a  result  of  the  random  dispersion  of  radial  orientations.  The  1 
axis  lies  along  the  fiber  axis,  but  the  23  plane  orientation  will  be  different  for  each  fibril. 
In  order  to  convert  the  raw  results  into  transverse  isotropic  material  properties,  the 
average  of  E2  and  E3  were  calculated  and  the  average  of  Gn  and  Gu  were  calculated. 

Additionally,  only  two  Poisson’s  ratios  are  needed;  these  were  found  by  averaging  v21 

and  v31  to  find  v2l ,  while  v32  and  v23  were  averaged  to  find  v32 .  The  resulting  values  are 
shown  in  Table  7. 


Table  7.  Transverse  isotropic  fibril  results.  Young’s  modulus  and  shear 

modulus  shown  in  GPa. 


Model 

e ; 

F 

-^23 

G23 

^12 

^21 

^32 

Linear 

Compression 

6.050 

6.601 

1.534 

1.677 

0.299 

0.295 

Tension 

5.760 

6.298 

1.458 

1.593 

0.299 

0.295 

Twisting 

Compression 

4.264 

3.812 

1.354 

1.368 

0.309 

0.366 

Tension 

4.054 

3.621 

1.286 

1.299 

0.309 

0.366 

These  results  can  be  compared  to  the  experimental  and  theoretical  values  shown 
in  Table  8.  The  experimental  results  of  nanoindentation  and  dynamic  mechanical  analysis 
evaluate  the  transverse  stiffness  of  a  fibril  in  compression.  The  other  results,  unless  noted, 
tested  the  longitudinal  stiffness  in  tension. 
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Table  8.  Comparison  of  values  for  fibril  results. 


Stiffness 

[GPa] 

Method 

Source 

7.65±3.85 

Nanoindentation  of  rat  tail  tendon  (Transverse) 

[56] 

1.00±0.75 

Finite  element  model,  tension  (Longitudinal) 

[11] 

7.5±5.5 

Finite  element  model,  compression  (Longitudinal) 

[11] 

4.96±0.57 

Dynamic  mechanical  analysis  using  AFM  ,  peak 
(Transverse)1 

[58] 

3.07±0.23 

Dynamic  mechanical  analysis  using  AFM,  hough 
(Transverse)1 

[58] 

4.36 

Molecular  multi-scale  modeling  (small  strain) 

[29] 

38 

Molecular  multi-scale  modehng  (large  strain) 

[29] 

10±2.0 

3D  molecular  dynamics  (Longitudinal) 

[48] 

5 

Molecular  modeling  (Longitudinal) 

[50] 

4.81 

Molecular  dynamics  and  finite  element  modehng 
(Longitudinal) 

[54] 

2.05±0.75 

Full  atomistic  model  (Longitudinal) 

[19] 

4.75±3.06 

Finite  element  model  (Longitudinal) 

[55] 

'Peak  and  trough  refer  to  mineralized  and  un-mineralized  sections,  respectively. 


As  the  results  shown  in  Table  8  confirm,  the  micromechanics  model  utilized  for 
bone  fibrils  is  within  the  range  of  values  defined  for  bone  fibrils.  The  model  provides  a 
simple  solution  to  a  complex  problem.  The  relatively  elementary  approach  yielded 
realistic  results  for  both  the  linear  fibril  model  and  the  twisting  fibril  model.  Additionally, 
the  micromechanics  model  calculated  the  Poisson’s  ratios  for  the  fibril  rather  than 
assuming  a  value,  as  did  all  results  in  Table  8. 

4.  Fibril  Discussion 

The  micromechanics  model  produces  accurate  results  for  the  stiffness  of  bone 
fibrils  in  tension  and  compression.  The  calculation  of  transverse  isotropic  material 
properties,  to  include  Poisson’s  ratios  and  shear  moduli,  can  be  used  to  increase  the 
accuracy  of  subsequent  hierarchies,  as  the  micromeclianics  model  will  be  used  once  again 
for  the  calculation  of  fiber  properties. 

One  metric  for  comparing  validity  of  theoretical  models  is  to  compare  bone 
mineral  content.  However,  veiy  little  information  is  available  as  to  the  mineral  content  of 
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microstructures.  Therefore,  the  mineral  content  at  each  level  will  be  tracked  for 
comparison  at  the  macro  level.  Table  9  displays  the  mineral  content  of  the  two  fibril 
models  proposed. 


Table  9.  Fibril  model  mineral  volume  fraction. 


Model 

Mineral  Volume  Fraction  [%] 

Linear 

16.66 

Twisting 

6.22 

B.  BONE  FIBER 

The  successive  hierarchy  following  the  bone  fibril  is  that  of  the  bone  fiber.  Fibers 
are  a  composite  material  of  fibrils  and  HA.  The  HA  is  present  in  the  fibers  in  the  form  of 
extra  fibrillar  mineral,  deposited  between  densely  packed  fibrils  [2],  [9],  [14],  [62], 

1.  Fiber  Structure 

Bone  fibers  are  sometimes  referred  to  as  fibril  bundles.  After  fibrils  have  been 
arranged  in  their  three-dimensional  space,  mineral  is  deposited  on  the  outside  of  the 
fibril.  This  extra  fibrillar  mineral  acts  as  a  stiffener  for  the  biurdles.  These  mineralized 
fibrils  are  tightly  packed  in  a  uniform  direction  [9],  similar  to  a  fiber  reinforced 
composite  [7],  This  grouping  of  fibrils  creates  the  fiber  structure.  The  mineralization  of 
fibrils  is  shown  in  Figure  1 1 . 
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Figure  1 1 .  Scanning  electron  micrograph  of  a  mineralized  fibril,  from  [63]. 

a.  Three-Dimensional  Fiber  Structure 

Each  closely  packed  fibril  is  surrounded  by  a  crust  of  mineral,  approximately  20- 
30  nni  thick  [7],  [62],  The  HA  is  grown  as  small  plates  that  combine  to  coat  the  fibrils 
[7],  [9].  This  provides  an  ordered  arrangement  of  minerals,  where  the  1  axis  of  minerals 
and  fibrils  run  parallel  to  each  other  [19],  [62],  [63].  This  allows  us  to  represent  HA  as  a 
transversely  isotropic  material.  The  thickness  of  the  crust  will  be  estimated  as  26  nm. 
Additionally,  a  single  uniform  fibril  diameter  of  150  nm  will  be  assumed  [7],  [9],  This 
three-dimensional  structure  is  shown  in  Figure  12,  with  the  surface  mineral  removed. 
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Figure  12.  Three-dimensional  fiber  array.  Dark  rectangles  represent  FIA  crystals 

while  the  rod  structures  represent  fibrils. 

The  fibrils  of  Figure  12  exhibit  a  linear  crystal  pattern,  however,  both  linear  and 
twisting  mineral  patterns  will  be  analyzed.  Very  little  is  known  about  the  degree  of 
mineralization  along  the  lengths  of  the  fibrils  [3],  [9].  In  order  to  analyze  the  effects  of 
extrafibrillar  mineralization  (EFM),  several  degrees  of  mineralization  will  be  tested. 

b.  Inclusion  of  N on-structural  Materials 

Mineral  packing  at  the  fiber  level  is  not  a  close-packed  assembly,  as  Figure  12 
shows.  The  minerals  allow  for  void  spaces  between  them.  These  voids  at  the  micro  level 
allow  for  the  presence  of  water  with  dissolved  non-stractural  proteins  and 
macromolecules  [9].  This  solution  is  often  referred  to  as  the  extrafibrillar  matrix.  This 
liquid  phase  at  the  fiber  level  represents  the  remaining  component  of  bone  water  [11] 
These  spaces  will  be  treated  as  water,  resulting  in  a  bi-modulus  property  for  fibers. 

2.  Micromechanical  Fiber  Model 

The  fibers  found  in  the  macro  level  are  single  bundles  of  fibrils  or  a  combination 
of  many  bundles.  The  bundles  are  an  organized  and  repeating  arrangement  of  fibrils  and 
HA.  Accordingly,  the  micromechanics  model  is  aptly  suited  for  modeling  bone  fiber.  The 
fiber  model  will  not  rely  on  exact  dimensions,  as  did  the  unit  cell  for  the  fibril  models. 
The  fiber  unit  cells  will  rely  on  volume  fractions  of  each  material.  In  order  to  ensure  a 
volume  of  unity,  Equation  18  will  be  adhered  to. 
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(18) 


°l  +  a2  ~  b\  +  b2  =  Cj  +  c2  =  1 

Additionally,  since  all  fibrils  are  assumed  circular,  symmetry  can  be  applied.  This 
allows  a  mineralized  fibril  to  be  split  into  four  equal  quadrants. 

a.  Unit  Cell  Dimensions 

For  all  models,  the  longitudinal  cross  section  will  remain  constant,  as  derived 
from  the  area  percent  HA  on  the  1  lace.  This  area  percent  assumes  a  total  crust  thickness 
of  26  mn  and  a  fibril  diameter  of  150  mn.  When  symmetry  is  applied,  a  crust  of  13  lmi 
and  a  radius  of  75  urn  are  utilized.  In  order  to  apply  this  to  the  micromechanics  model, 
the  dimensions  shown  hi  Table  10  will  be  applied  to  the  1  face. 


Table  10.  Fiber  unit  cell  1  face  dimension. 


Dimension 

Value 

by 

85 

b2 

15 

cl 

85 

C2 

15 

This  equates  to  a  fibril  area  of  72.25%  when  converted  into  the  micromechanics 
model.  This  is  visualized  in  Figure  13. 
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Figure  13.  Longitudinal  fiber  cross  section  as  compared  to  micromechanics 
model  cross  section,  (a)  Fiber  cross  section:  dark  area  represents 
mineral,  light  section  represents  fibril,  (b)  Micromechanics  model 
longitudinal  cross  section. 


Subcell  1  represents  the  fibril.  In  order  to  test  the  effect  of  fibril  mineralization, 
the  degree  of  mineralization  will  be  varied.  %EFM  will  be  defined  as  the  total  percent  of 
fibril  covered  by  mineral  and  is  equal  to  the  volume  of  subcells  3,  5,  and  7,  shown  in 
Figure  14. 
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Figure  14. 


By  altering  the  length  of<7j ,  %EFM  is  modified  for  the  fiber.  We  will  assume  that 

the  fiber  is  comprised  of  no  more  than  95%  EFM.  The  different  levels  of  mineralization 
analyzed  are  listed  in  Table  1 1 . 


Table  1 1 .  Fiber  unit  cell  extrafibrillar  mineralization  dimensions. 


%EFM 

a2 

50 

50 

50 

70 

70 

30 

90 

90 

10 

95 

95 

5 
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In  addition  to  altering  the  %EFM  of  the  unit  cell,  the  effects  of  water  will  be  taken 
into  account.  In  compression,  subcells  4,  6,  and  8  will  be  assigned  the  properties  of  water. 
In  tension,  these  subcells  will  be  assigned  as  void  spaces.  In  order  to  apply  stiffness  and 
Poisson’s  ratio  to  water,  it  was  assumed  to  be  slightly  compressible.  Additionally,  the 
elastic  modulus  was  assiuned  equal  to  the  bulk  modulus  when  confined  to  small 
constrained  volumes,  as  with  microscale  bone.  Therefore,  the  Poisson’s  ratio  can  be 
solved  using  Equation  19. 

3  K-E  3  K-K  2  K  1  /in, 

v  = - = - = - =  -  (19) 

6K  6K  6K  3 

where  K  is  the  bulk  modulus.  The  material  properties  for  water  and  void  space  are  listed 
in  Table  12.  The  material  properties  for  HA  and  the  fibril  are  listed  in  Table  4  and  Table 
7,  respectively. 


Table  12.  Material  properties  of  water  and  void  space. 


Material 

E  [GPa] 

V 

G  [GPa] 

Water 

2.15 

0.33 

0.001 

Void 

0.001 

0.001 

0.001 

b.  Results 

The  results  for  the  fiber  model  were  calculated  both  in  compression  and  tension. 
As  expected,  due  to  the  applied  symmetry,  the  fiber  exhibits  transverse  isotropic  material 
properties.  This  is  due  to  the  assumption  that  the  radial  arrangement  of  the  fiber  arrays  is 
a  random  occurrence.  The  results  from  the  fiber  models  are  split  into  two  tables  for  easier 
comparison.  Table  13  shows  the  results  in  compression,  while  Table  14  shows  the  results 
in  tension. 
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T able  13.  Fiber  results  in  compression. 


50°/ 

Linear 

o  EFM 
Twisting 

70°/ 

Linear 

>EFM 

Twisting 

90°/ 

Linear 

o  EFM 
Twisting 

95°/ 

Linear 

o  EFM 
Twisting 

El  [GPa] 

6.511 

5.446 

7.567 

6.509 

11.60 

10.47 

16.17 

14.97 

E23  [GPa] 

16.66 

14.62 

21.41 

19.16 

26.35 

23.91 

27.72 

25.23 

G13  [GPa] 

1.053 

0.931 

1.473 

1.302 

1.893 

1.672 

1.998 

1.765 

Gn  [GPa] 

0.012 

0.012 

0.019 

0.019 

0.056 

0.056 

0.109 

0.108 

*21 

0.369 

0.413 

0.364 

0.397 

0.295 

0.307 

0.245 

0.250 

*32 

0.140 

0.132 

0.128 

0.119 

0.126 

0.118 

0.127 

0.120 

The  results  in  compression  show  that  the  linear  fibril  model  produces  a  stiffer 
fiber  in  all  normal  and  shear  cases  and  in  all  directions.  For  all  %EFM  the  fiber  shear 
moditlus  in  the  12  direction  is  much  less  than  the  23  direction.  As  %EFM  increases  from 
50%  to  95%,  the  12  shear  modulus  increases  by  an  order  of  magnitude,  while  the  23 
modulus  approximately  doubles.  This  is  shown  for  both  the  linear  and  the  twisting 
models.  The  Poisson’s  ratio  in  the  32  dir  ection  remains  relatively  constant,  but  as  %EFM 
increases  the  Poisson’s  ratio  in  the  21  direction  decreases. 


Table  14.  Fiber  results  in  tension. 


50°/ 

Linear 

o  EFM 
Twisting 

70°/( 

Linear 

>EFM 

Twisting 

90°/ 

Linear 

EFM 

Twisting 

95°/ 

Linear 

o  EFM 
Twisting 

El  [GPa] 

4.527 

3.374 

4.691 

3.592 

4.870 

3.841 

4.919 

3.911 

-^23  [GPa] 

14.10 

12.80 

19.74 

17.95 

25.38 

23.03 

26.79 

24.31 

G23  [GPa] 

1.001 

0.885 

1.401 

1.237 

1.800 

1.589 

1.900 

1.677 

0,2  [GPa] 

0.012 

0.012 

0.019 

0.019 

0.056 

0.056 

0.109 

0.107 

*21 

0.227 

0.266 

0.318 

0.372 

0.408 

0.478 

0.431 

0.505 

*32 

0.113 

0.099 

0.113 

0.099 

0.113 

0.099 

0.113 

0.099 

The  results  in  tension  show  similar  results  to  those  in  compression.  The  linear 
model  exhibits  greater  stiffness  for  all  results.  Additionally,  the  values  of  E23  are  much 
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greater  than  Ex .  The  value  of  v32  remains  constant  for  each  model  as  %EFM  increases. 
However,  the  value  of  v2i  increases  steadily  as  %EFM  increases. 

3.  Fiber  Results  Discussion 

The  shear  moduli  in  the  12  direction  are  very  low  for  small  %EFM.  This  is  due  to 
the  void  relatively  large  amount  of  void  space  present  within  the  unit  cell.  However,  this 
will  be  redressed  with  structure  of  lamellar  bone.  The  inclusion  of  fibers  within  the 
disordered  fibrillar  matrix  will  increase  the  shear  modulus  in  the  12  direction. 

The  difficulty  with  evaluating  the  results  of  the  fiber  model  is  that  there  is  an 
absence  of  experimental  testing  available  for  comparison.  Fibers  exist  within  the 
macrostructures  of  bone  and  are  not  easily  isolated  for  testing.  Additionally,  almost  all 
theoretical  calculations  assume  that  macrostructure  bone  is  composed  solely  of  layered 
fibers.  Due  to  the  recent  findings  of  [12],  we  know  this  not  to  be  true.  The 
macrostructures  of  bone  are  a  fiber  reinforced  composite  with  fibrils  acting  as  the  matrix 
and  bone  fibers  as  the  fibers. 

The  mineral  content  completely  surrounds  the  fibril.  This  complete  encirclement 
increases  the  normal  stiffness  in  the  23  direction  and  stiffens  the  fiber  against  shear-  on 
the  23  plane.  The  mineral  content  of  the  different  fiber  models  was  calculated.  These 
calculations  include  both  intra fibrillar  and  extrafibrillar  HA.  The  resulting  vahres  are 
shown  in  Table  15. 


Table  15.  Fiber  model  mineral  volume  fraction. 


%EFM 

Mi  tier  a  1  \ 
Fraction 
Linear 

olume 

[%] 

Twisting 

50 

25.91 

18.37 

70 

31.46 

23.92 

90 

37.01 

29.47 

95 

38.40 

30.86 
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C.  MICROSCALE  CONSIDERATIONS 

The  models  presented  for  the  microscale  components  of  bone  provide  simple  and 
accurate  solutions.  They  calculate  transverse  isotropic  material  properties  driven  by  the 
ordered  geometry  inherent  to  the  self-assembly  of  living  materials.  However,  not 
everything  is  known  about  the  microstructures  of  bone. 

1.  Biological  Unknowns 

There  is  very  little  information  as  to  the  effective  length  of  microstructures.  It  is 
difficult  to  produce  an  image  detailing  the  beginning  and  end  of  fibrils  or  fibers  [62],  In 
addition  to  not  knowing  the  length  of  these  structures,  little  information  is  known  about 
the  attachment  of  these  endpoints  or  how  they  are  secured  within  the  surrounding 
material.  The  concept  of  twisting  HA  crystals  in  fibrils  is  a  recent  proposition  [9],  As 
Figure  5b  shows,  the  lateral  spacing  is  not  restricted  to  a  67  nm  periodicity.  With  no 
definitive  model  of  three-dimensional  fibril  packing,  the  linear  fibril  model  presents  the 
most  biologically  probable  result. 

The  micromechanics  model  developed  relies  on  selecting  the  proper  unit  cell 
dimensions.  With  the  linear  model,  the  dimensions  presented  are  logically  derived.  The 
twisting  model  is  based  on  the  theoretical  assumption  that  an  incorporated  twist  is  a 
natural  occurrence.  This  introduces  two  unknowns  that  need  to  be  addressed.  The  rate  of 
twist  is  undefined  and  the  twist  is  not  restricted  to  lateral  packing,  but  to  the  surface  of  a 
fibril  with  a  circular  cross  section.  Thus,  the  models  presented  in  this  study  assumed  a 
lateral  packing  with  a  periodicity  of  67nm.  This  introduces  a  twist  in  the  HA  crystal 
pattern.  Due  to  these  assumptions  the  unit  cell  dimensions  are  less  quantitatively  defined 
for  the  twisting  fibril  model. 

Furthermore,  the  effect  of  shearing  between  fibrils  and  fibers  is  a  currently 
researched  topic.  Several  finite  element  models  and  molecular  dynamics  models  have 
attempted  to  explain  these  occurrences  through  the  use  of  bonding  between  different 
phases  [7],  [54]  and  energy  dissipation  of  cross-links  [29],  [50].  However,  no  agreed 
upon  mechanisms  have  been  produced  through  these  computationally  expensive  finite 
element  methods  (FEM)  and  molecular  dynamics  models. 
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These  models  hypothesize  the  chemical  nature  of  in  vivo  bone  at  the  nano  and 
micro  level.  This  information  is  difficult  to  obtain  through  laboratory  testing.  The 
challenge  is  understandable,  as  recreating  the  in  vivo  conditions  in  a  laboratory  are 
extremely  difficult.  In  addition,  the  living  cells  present  in  micro  structures  provide 
continuous  growth  and  remodeling  of  bone  structures.  These  living  cells  provide  another 
unknown  variable  to  the  modeling  of  bone  structures. 

Where  information  is  lacking,  assumptions  allow  results  to  be  obtained.  The 
models  presented  here  help  to  make  simple  assumptions  that  allow  for  the  lack  of  detailed 
information  to  be  negated.  As  a  result,  the  models  match  current  hypotheses  on  the 
structure  of  microscale  bone  components. 

2.  Mineral  Content 

An  important  metric  for  evaluating  bone  is  the  mineral  content.  Studies  have  been 
completed  comparing  volume  percent  and  weight  percent  of  HA  present  in  different 
structures  of  bone.  However,  many  of  these  studies  are  at  the  macro  level.  A  comparison 
will  be  made,  once  the  macro  level  models  are  complete,  of  the  various  levels  of 
mineralization. 


3.  Variations  in  Experimental  Testing 

As  Table  8  shows,  the  experimental  and  theoretical  results  for  biological  materials 
are  widely  varied.  Variations  in  theoretical  results  are  due  to  the  assumptions  made  at 
each  level,  the  scope  of  the  testing,  the  source  of  material  properties,  the  complexity  of 
each  model,  etc.  Variations  in  experimental  results  are  due  to  the  conditions  of 
experimental  tests  and  the  methods  used  for  each  experimental  test.  Experimental  results 
test  hydrated  structures,  dehydrated  structures,  structures  in  vitro,  in  vivo  structures 
prepared  in  a  series  of  solutions,  non-human  bone  structures  and  many  other  variations. 
With  so  many  variations,  how  is  it  possible  to  select  the  most  accurate  result? 

It  is  for  this  reason  that  the  linear  and  twisting  fibrillar  models  were  both  used  in 
the  fiber  model.  It  is  the  same  reason  that  multiple  fiber  models  will  be  applied  to  the 
macro  level  models.  These  results  can  then  be  compared  to  experimental  results  of  human 
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bone.  The  results  of  macro  scale  bone  are  more  easily  standardized  and  tested,  leading  to 
more  uniform  experimental  results.  Additionally,  once  the  macro  scale  models  are 
compared  to  experimental  results,  the  accuracy  of  each  microscale  model  can  be 
retroactively  verified. 
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IV.  MACROSCALE 


From  the  basic  building  blocks  of  HA  and  TC,  complex  microscale  geometries 
have  been  created.  These  microscale  components  serve  as  the  elementary  materials  for 
further  biocomposite  materials.  In  order  to  analyze  the  macroscale  structures  of  bone, 
three  different  models  will  be  utilized:  a  micromechanics  model  of  lamellar  bone,  a 
layered  fiber  reinforced  composite  model  of  cortical  and  trabecular  bone,  and  a  finite 
element  model  of  a  tetrakaidecahedron  modeling  cancellous  bone. 

A.  LAMELLAR  BONE 

The  next  hierarchical  step  of  bone  is  that  of  lamellar  bone.  Lamellar  bone  is  a 
fiber  reinforced  composite  of  bone  fibers  and  bone  fibrils.  The  presence  of  a  disordered 
fibril  matrix  is  a  relatively  new  discovery,  identified  through  the  use  of  a  serial  surface 
view  (SSV).  By  continuously  etching  and  imaging  a  section  of  material,  the  three 
dimensional  structure  is  revealed.  Reznikov  et  al.  [12]  has  utilized  SSV  to  three 
dimensionally  deconstruct  sections  of  bone. 

1.  Lamellar  Structure 

Lamellar  bone  is  comprised  of  bone  fibers  surrounded  by  a  matrix  of  disordered 
fibrils  [9],  [12],  This  is  a  new  revelation  that  many  older  models  do  not  take  into  account. 
Previous  calculations  have  assumed  that  lamellar  layers  are  merely  an  ordered  array  of 
fibers  [1],  [3],  [10],  [13],  [64],  The  discovery  of  the  disordered  fibril  matrix  adds  another 
hierarchical  level  that  can  be  modeled  using  the  micromechanics  model.  Reznikov  et  al. 
[12]  found  that  the  lamellar  bone  exhibits  an  ordered  motif  of  2-3  jum,  representing  the 
layer  thickness,  and  disordered  matrix  thickness  of  0.25  /um  to  1  /um,  representing  the 
fibrillar  matrix.  In  addition  to  the  presence  of  a  disordered  matrix,  the  SSV  discovered  the 
presence  of  small  pores  with  a  diameter  of  50  nm  [12],  These  pores  are  represented  as 
voids  in  the  lamellar  matrix,  to  be  ignored  in  this  model.  The  simplified  view  of  a  single 
lamellar  layer  is  shown  in  Figure  15. 
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1 

Figure  15.  Single  lamellar  bone  layer. 

Although  still  measured  in  microns,  the  lamellar  layer  is  quantified  in  the 
macroscale  section  due  to  the  fact  that  it  is  the  main  constituent  of  cortical  and  cancellous 
bone. 


2.  Lamellar  Layer  Micromechanics  Model 

The  calculation  of  material  properties  for  lamellar  bone  will  closely  mirror  the 
process  used  for  analyzing  the  material  properties  of  bone  fibers.  Utilizing  the  knowledge 
that  fibers  within  a  lamellar  layer  are  unidirectional  and  that  the  disordered  matrix  has  a 
relatively  random  thickness,  within  the  bounds  defined  by  [12],  symmetry  can  be  applied 
to  solve  for  transverse  isotropic  material  properties.  The  subunit  to  be  used  for  the 
lamellar  layers  will  assume  that  the  layer  thickness  is  2.5  /um  and  the  matrix  thickness  is 
0.375  fim.  These  are  the  average  of  the  range  defined  by  [12],  The  fiber  diameter  is  then 
calculated  to  be  1 .75  /um .  The  longitudinal  cross  section  to  be  used  for  lamellar 
calculations  is  shown  in  Figure  16. 
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(a) 


Figure  16.  Longitudinal  lamellar  cross  section  as  compared  to  micromechanics 
model  cross  section,  (a)  Lamellar  cross  section:  shaded  area  represents 
disordered  matrix,  light  section  represents  fiber,  (b)  Micromechanics 
model  longitudinal  cross  section. 


The  longitudinal  cross  section  of  the  micromechanics  model,  shown  in  Figure 
16b,  provides  a  38.4%  fiber  volume  for  the  lamellar  unit  cell.  This  is  derived  from  the 
assumed  fiber  diameter  and  matrix  thickness  between  fibers.  The  presence  of  the  small 
pores  will  be  neglected  as  they  account  for  approximately  0.1%  of  the  current  model’s 
volume  and  are  not  uniformly  present  along  the  length  of  the  fibers. 

The  micromechanics  model  will  once  again  utilize  volume  percent  to  analyze  the 
lamellar  layer;  unit  cell  dimensions  will  be  constrained  to  Equation  18.  In  order  to  model 
the  lamellar  layer  as  a  continuous  fiber,  the  unit  cell  shown  in  Figure  17  is  used. 
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Figure  17.  Unit  cell  for  lamellar  micromechanics  model. 


The  material  properties  assigned  to  subcells  1,  3,  5  and  7  are  mirrored  for  subcells 
2,  4,  6  and  8.  This  produces  a  continuous  fiber  reinforced  composite. 

a.  Unit  Cell  Quantification 

The  dimensions  of  ax  and  a2  are  unimportant,  as  the  properties  of  the  odd 

subcells  are  mirrored  by  the  even  subcells.  The  dimensions  to  be  used  for  this  model  are 
listed  in  Table  16. 
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Table  16.  Unit  cell  dimension  for  lamellar  model. 


Dimension 

Value 

ai 

50 

ai 

50 

62 

h2 

38 

CX 

62 

C2 

38 

The  composite  matrix  is  composed  of  a  disordered  array  of  fibrils.  The  fibrils  do 
not  exhibit  any  preferential  alignment  [12],  therefore  the  material  can  be  assumed 
isotropic.  This  will  be  accomplished  by  averaging  the  stiffnesses  and  Poisson’s  ratios 
found  in  Table  7.  Furthermore,  fiber  models  of  varying  %EFM  will  be  used.  Due  to  the 
bi-modulus  property  of  the  materials  within  the  model  and  the  applied  symmetry  of  the 
unit  cell,  the  lamellar  layer  is  expected  to  show  transverse  isotropic,  bi -modulus,  material 
properties. 

b.  Results 

The  results  of  the  lamellar  ruicromechanics  model  were  calculated  in  both 
compression  and  tension  for  four  different  %EFM.  The  resulting  stiffnesses  and 
Poisson’s  ratios  are  shown  in  Table  17  and  Table  18. 
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Table  17.  Lamellar  results  in  compression. 


50°/< 

Linear 

EFM 

Twisting 

70°/ 

Linear 

EFM 

Twisting 

90°A 

Linear 

>EFM 

Twisting 

95°/ 

Linear 

o  EFM 
Twisting 

Ex  [GPa] 

6.428 

4.609 

6.841 

5.024 

8.391 

6.547 

10.14 

8.272 

F 

-'-'23 

[GPa] 

9.001 

6.425 

9.717 

6.928 

10.35 

7.409 

10.55 

7.576 

G23 

[GPa] 

1.336 

1.156 

1.552 

1.338 

1.705 

1.466 

1.737 

1.492 

<*u 

[GPa] 

0.622 

0.529 

0.629 

0.536 

0.665 

0.572 

0.715 

0.620 

V2l 

0.328 

0.366 

0.322 

0.350 

0.280 

0.286 

0.242 

0.237 

v32 

0.237 

0.266 

0.235 

0.269 

0.243 

0.288 

0.252 

0.302 

The  results  of  Table  17  show  that  the  linear  model  exhibits  greater  stiffness  than 
the  twisting  model  for  all  instances.  Additionally,  for  low  %EFM  the  Yoimg’s  modulus 
in  the  1  direction  is  smaller  than  in  the  23  direction.  As  %EFM  increases,  the  two  values 
converge.  For  the  twisting  model,  the  stiffness  in  the  1  direction  surpasses  the  stiffness  in 
the  23  direction. 


T able  1 8 .  Lamellar  results  in  tension. 


50°/< 

Linear 

>EFM 

Twisting 

70°/i 

Linear 

>EFM 

Twisting 

90°/ 

Linear 

.EFM 

Twisting 

95°/ 

Linear 

0  EFM 
Twisting 

El  [GPa] 

5.514 

3.718 

5.578 

3.801 

5.646 

3.895 

5.665 

3.922 

F 

£,23 

[GPa] 

8.305 

5.992 

9.174 

6.538 

9.781 

6.902 

9.905 

6.976 

G23  [GPa] 

1.270 

1.098 

1.475 

1.271 

1.621 

1.393 

1.651 

1.418 

Gn 

[GPa] 

0.592 

0.503 

0.599 

0.510 

0.635 

0.546 

0.684 

0.593 

V2\ 

0.312 

0.370 

0.338 

0.394 

0.355 

0.407 

0.359 

0.410 

v32 

0.235 

0.262 

0.230 

0.257 

0.226 

0.254 

0.225 

0.253 
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The  results  in  tension  also  show  that  the  linear  model  provides  greater  stiffness 
than  all  twisting  models.  Additionally,  the  stiffness  in  the  23  direction  is  larger  than  the 
stiffness  in  the  1  dir  ection  for  all  %EFM. 

3.  Lamellar  Discussion 

As  with  microscale  hierarchies,  the  mineral  content  of  the  lamellar  layers  was 
analyzed.  The  mineral  volume  fraction  of  each  model  is  displayed  in  Table  19. 


Table  19.  Lamellar  model  mineral  volume  fraction. 


%EFM 

Mineral  Vo 
Linear 

ume  Fraction  [%] 
Twisting 

50 

20.22 

10.89 

70 

22.35 

13.02 

90 

24.48 

15.16 

95 

25.02 

15.69 

The  mineral  volume  fraction  within  lamellar  layers  is  representative  of 
macroscale  bone  mineral  content.  Comparisons  can  be  made  to  theoretical  and 
experimental  values.  Early  studies  have  hypothesized  on  the  mineral  content  of  bone 
through  both  estimations  of  spatial  arrangement  and  calculation  from  bone  ash  content. 
Early  understandings  of  fibril  organization  place  total  bone  mineral  volume  content  at 
50%  [65].  Calculation  of  the  bone  ash  content  of  adult  cows  was  found  to  be  close  to 
70%  [66],  The  volume  percent  of  bone  ash  does  not  correlate  to  bone  mineral  content. 
Bone  ash  is  determined  by  burning  a  sample  of  bone  in  a  furnace  and  comparing  the  ash 
weight.  This  method  results  in  a  large  estimate  of  mineral  content  as  additional  residues 
remain  from  sources  other  than  hydroxyapatite.  More  recent  studies  estimate  a  lower 
bone  mineral  volume.  Kotha  and  Guzelsu  [67]  postulate  that  mineral  volume  content  is 
40%  for  bone  while  [68]  approximates  33-43%  apatite  mineral  volume.  These  estimates 
are  based  on  a  more  complete  understanding  of  the  hierarchical  structure  of  bone,  but  are 
still  estimates  as  bone  mineral  content  varies  with  bone  type,  anatomical  position,  age, 
and  gender. 
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Comparing  the  bone  mineral  content  results  of  this  study  to  those  found  through 
experimental  and  theoretical  methods,  the  linear  crystal  pattern  emerges  as  the  more 
viable  model.  The  mineral  volume  fractions  of  the  twisting  model  are  too  low  to  validate 
its  use.  Even  with  slight  perturbations  to  the  constraints  applied  to  the  fibril  and  fiber 
models,  the  mineral  content  of  the  twisting  model  does  not  match  the  current  estimates. 
Future  calculations  will  assume  the  linear  model  to  be  the  valid  model. 

B.  CORTICAL  BONE 

Cortical  bone  is  the  dense  structure  of  macroscale  bone  [6],  It  provides  an  external 
shell  surrounding  the  spongy  center  of  cancellous  bone  and  accounts  for  a  large  percent 
of  the  weight  fraction  of  the  skeletal  system.  Cortical  bone  is  composed  of  concentric 
layers  of  lamellar  bone,  known  as  osteons.  Osteons  can  be  made  of  primary  or  secondary 
bone.  Primary  bone  is  found  where  bone  has  been  grown  de  novo.  Secondary  bone  is 
created  when  the  primary  bone  is  broken  down  and  re-grown  in  place  [9].  Secondary 
osteons  are  called  Haversian  systems,  which  will  be  analyzed  in  this  study. 

1.  Three  Dimensional  Haversian  System  Structure 

The  concentric  layers  composing  Haversian  systems  form  a  cylindrical  structure 
approximately  200  microns  in  diameter  [9],  [12],  [64],  [69].  At  the  center  of  these 
cylinders  is  a  Haversian  canal.  These  canals  contain  the  blood  supply  and  nerve  endings 
for  the  surrounding  bone.  The  canals  are  approximately  30-50  microns  in  diameter  [9], 
[12],  [64],  [69], 

Each  concentric  layer  associated  with  the  osteon  is  composed  of  lamellar  bone. 
The  fibers  of  each  individual  layer  are  unidirectional  in  alignment  [9],  [12],  [64],  [69]. 
The  center  cylinder  is  representative  of  a  Haversian  canal  and  will  be  represented  as  un¬ 
bound  water  [9].  Multiple  Haversian  systems  are  packed  together  in  cortical  bone.  Due  to 
their  circular  shape,  there  are  incomplete  layers  at  the  interface  of  each  Haversian  system 
where  the  boundaries  intersect.  These  boundaries  are  defined  by  a  cement  line,  which  is 
an  identifiable  region  where  osteon  growth  direction  has  transitioned.  The  properties  of 
cement  lines  are  similar  to  those  of  the  surrounding  bone,  despite  the  misnomer  of 
“cement”  [70], 
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The  long  axis  of  a  Haversian  system  is  preferentially  aligned  with  the  long  axis  of 
the  bone  for  a  majority  of  cortical  bone  [64],  Transversely  oriented  Haversian  canals  do 
exist  that  interconnect  each  longitudinal  Haversian  canal,  but  they  are  very  short 
compared  to  the  longitudinal  canals.  This  study  will  assume  all  Haversian  systems  are 
aligned  parallel  to  the  long  axis  of  the  bone.  Additionally,  as  the  osteoclasts  and 
osteoblasts  perform  the  resorption  and  deposition  process  they  leave  behind  a  small 
network  of  voids.  These  voids  are  known  as  canaliculi  and  promote  the  metabolic 
function  of  bone  by  transporting  proteins  and  molecules  within  the  densely  packed 
osteons  [71]. 

The  fiber  directions  of  the  concentric  lamellar  layers  vary  for  each  layer.  Two 
distinct  patterns  have  been  identified  in  lamellar  bone:  a  periodic  alternating  pattern  and  a 
continuous  fiber  twist.  The  alternating  pattern  was  found  by  both  [12]  and  [64],  Reznikov 
et  al.  [12]  utilize  SSV,  while  Varge  et  al.  [64]  used  synchrotron  X-ray  phase  nano- 
tomography  (SR-PNT).  SR-PNT  works  by  taking  different  tomographic  images  at 
periodic  distances  within  a  sample.  Both  sources  identified  a  distinct  alternating  pattern 
of  high  and  low  angle  fibers.  This  motif  exhibits  a  layer  of  fibers  at  65-80°  and  a  layer  at 
15-30°.  Reznikov  et  al.  [12]  refers  to  these  as  transverse  and  longitudinal  orientations, 
respectively.  Additionally,  due  to  the  imperfect  nature  of  biomaterials,  an  additional 
percent  of  intermediate  angles  were  found.  The  average  orientations  of  three  samples 
were  calculated.  The  samples  were  taken  from  a  20  year-old  female,  a  59  year-old  male, 
and  a  77  year-old  male.  The  results  revealed  45%  transverse  fibers,  35%  longitudinal 
fibers  and  20%  intermediate  fibers  [12],  The  findings  of  [64]  also  describe  a  section  of 
femur  exhibiting  a  continuous  twist  of  fibers  from  0-180°  between  lamellar  layers.  This 
continuous  twist  changes  by  approximately  10°  each  layer.  Both  the  periodic  alternating 
pattern  and  the  continuous  twist  can  be  modeled  as  a  composite  of  laminated  transverse 
isotropic  plates. 
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2.  Layered  Composite  Model 

The  layered  composite  of  macroscale  bone  is  composed  of  concentric  layers  of 
lamellar  bone.  The  previous  section  derived  the  properties  for  a  single  layer.  Two  layered 
motifs  will  be  analyzed  with  respect  to  cortical  bone. 

a.  Cortical  Bone  Model  Derivation 

The  stratified  lamellar  layers  of  cortical  bone  represent  a  combination  of 
laminated  plates.  Each  plate  is  represented  as  a  transverse  isotropic  material.  Due  to  the 
close  packing  of  the  lamellar  layers  within  Haversian  systems,  the  material  properties  of 
the  composite  can  be  found  by  analyzing  each  fiber  orientation.  The  breakdown  of  an 
osteon  into  its  respective  layers  can  be  seen  in  Figure  18. 


Figure  18.  Osteon  cross  sectional  view. 
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The  SSV  and  SR-PNT  views  of  cortical  bone  fiber  orientation  do  not  specify 
which  osteonal  layer  contains  what  fiber  direction  [12],  [64],  The  studies  are  a  focused 
view  of  a  single  section  of  bone.  Additionally,  due  to  the  close-packed  nature  of  adjacent 
Haversian  systems,  the  overall  material  properties  will  be  found  through  fractional 
volume  percent  of  each  fiber  orientation. 


The  material  properties  of  a  single  lamellar  layer  can  be  found  by  rotating  the 
layer  with  respect  to  the  3 -axis,  as  represented  in  Figure  15.  This  rotation  can  be 
completed  by  rotating  the  stiffness  matrix  of  the  material.  The  stiffness  matrix  is  a 
combined  view  of  the  generalized  Hooke’s  law.  The  two  derivations  of  Hooke’s  law  are 
shown  in  Equations  20-21  [72], 
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(20) 


(21) 


The  c  and  s  matrices  are  the  stiffness  and  compliance  matrices  of  the  material, 
respectively.  Due  to  the  symmetric  nature  of  the  lamellar  layers,  the  stiffness  and 
compliance  matrices  are  reduced  to  6  independent  constants.  The  stiffness  matrix  is 
shown  below. 
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C,1 

C12 

Cj2 

0 

0 

Cl  2 

C22 

C23 

0 

0 

Ci  2 

C23 

C22 

0 

0 

0 

0 

0 

C44 

0 

0 

0 

0 

0 

C66 

0 

0 

0 

0 

0 

(22) 


The  C„  and  Su  matrices  are  the  inverse  of  each  other.  In  order  to  enact  a  rotation 

of  each  lamellar  layer,  the  stiffness  matrix  will  be  rotated.  The  stiffness  matrix  can  be 
found  by  inverting  the  compliance  matrix  of  each  lamellar  layer.  The  normal  and  shear 
properties  were  derived  independently  for  each  material,  allowing  the  compliance  matrix 
to  be  built. 


In  order  to  solve  for  the  stiffness  matrix,  Equation  23  is  inverted.  This  matrix  can 
then  be  rotated  according  to  a  newly  defined  set  of  axis.  Equation  24  defines  the  rotated 
material  properties. 

o'  =  C'  s'  (24) 

The  prime  notation  indicates  the  newly  rotated  axis  coordinates,  where  C"  is  the 
stiffness  matrix  of  the  original  material.  In  order  to  solve  for  the  C  matrix  along  the 
original  axes,  the  stress  and  strain  must  be  transformed. 

7>  =  C'7>  (25) 
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By  inverting  Ta  and  pre-multiplying  into  the  right  side  of  Equation  25,  the 
generalized  Hooke’s  equation  is  found. 

<7  =  T~xC'Tes  (26) 

Equation  26  shows  that  the  stiffness  matrix  along  the  new  coordinates  can  be  defined  as: 

C  =  T-'C'T:  (27) 

where  T~] 1  is  the  inverse  of  the  stress  transformation  matrix  and  Te  is  the  transformation 

matrix  of  the  engineering  strain.  The  transformation  matrices  are  defined  in  Equations 
28-29. 

m 2  n2  0  0  0  -2 mn 

n2  m2  0  0  0  2 mn 

0  0  1  0  0  0 

0  0  0  in  n  0 

0  0  0  —n  m  0 

mn  — mn  0  0  0  ( /» 2  -  n 2  j 

m2  n2  0  0  0  mn 

n2  m2  0  0  0  —mn 

0  0  1  0  0  0 

0  0  0  m  —n  0 

0  0  0  n  m  0 

-2 mn  2 mn  0  0  0  ( m 2  -  n 2  j 

m  =  cos#  and  n  =  sin#,  where  6  is  defined  as  rotation  about  the  3  axis. 

In  addition  to  accounting  for  the  rotation  of  the  individual  lamellar  layers,  the 
Haversian  canal  and  canaliculi  will  be  addressed  in  the  macroscale  model.  Based  on  the 
dimensions  of  the  Haversian  system,  the  Haversian  canal  and  the  microscopic  canaliculi, 
macroscale  cortical  bone  is  75%  densely  packed  bone  and  25%  void  space.  As  with 
smaller  scale  models  the  void  space  will  be  addressed  as  a  liquid  in  compression  and  a 
void  in  tension. 
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b.  Fiber  Direction 

In  order  to  calculate  the  composite  properties  of  cortical  bone,  two  models  will  be 
assessed.  One  model  utilizes  the  preferential  orientations  found  by  [12]  while  the  other 
model  assumes  a  continuous  twist  as  found  by  [64],  These  two  models  will  be  referred  to 
as  preferential  orientation  and  smooth  orientation,  respectively,  from  this  point  on. 

The  two  models  rely  on  volume  percent  of  a  defined  fiber  orientation.  The  volume 
percentages  present  for  each  fiber  orientation  of  the  preferential  model  are  shown  in 
Table  20  and  the  designated  volume  percentages  for  the  smooth  model  are  shown  in 
Table  21. 


Table  20.  Preferential  orientation  layered  composite  model  parameters. 


Fiber  Designation 

Fiber  Dir  ection 

Volume  Percent 

Transverse 

72.5° 

45 

Longitudinal 

22.5° 

35 

Intermediate 

0° 

5 

Intermediate 

o 

O 

5 

Intermediate 

55° 

5 

Intermediate 

90° 

5 

The  transverse  and  longitudinal  fibers  are  set  at  the  average  of  the  range  defined 
by  [12],  The  intermediate  fibers  are  split  between  four  intermediate  directions,  at  5% 
volume  each. 


Table  2 1 .  Smooth  orientation  layered  composite  model  parameters. 


Fiber  Direction 

Volume  Percent 

0° 

10 

10° 

10 

20° 

10 

30° 

10 

O 

o 

10 

50° 

10 

60° 

10 

70° 

10 

80° 

10 

90° 

10 
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The  fiber  distribution  for  the  smooth  orientation  model  represents  a  consistent 
fanning  angle.  The  stiffness  matrix  of  each  layer  is  calculated  using  Equation  28,  using 
the  fiber  directions  listed  in  Table  20  and  Table  21  as  the  theta  values  for  Equations  29- 
30.  The  Voigt  estimation  of  composite  stiffness  is  calculated  for  the  volume  percent 
shown  in  Table  20  and  Table  21.  It  should  be  noted  that  the  rotation  of  fiber  dir  ections  is 
assumed  symmetrical  about  0°.  The  fiber  directions  as  stated  by  [12]  and  [64]  do  not 
designate  either  positive  or  negative  angles  about  the  longitudinal  axis.  Therefore,  the 
results  of  the  layered  composite  for  cortical  bone  model  are  expected  to  produce  an 
orthotropic  linear  elastic  material.  However,  due  to  the  random  radial  distribution  of 
osteons,  we  should  expect  a  transverse  isotropic  material.  Therefore,  the  average  of  the 
Sn  and  5,3  will  be  used  to  calculate  v21 ,  the  average  of  S22  and  S33  will  be  used  to  find 
E, 3  and  v32 ,  and  the  average  of  S55  and  S66  will  be  used  to  find  G12  . 

c.  Cortical  Model  Results 

The  results  of  the  layered  composite  model  were  calculated  using  only  the  linear 
fibril  model.  Additionally,  the  preferential  and  smooth  layered  models  were  calculated  in 
both  tension  and  compression.  The  compressive  results  are  shown  in  Table  22  while  the 
tensile  results  are  shown  in  Table  23. 


Table  22.  Cortical  bone  in  compression 


50% 

Pref. 

EFM 

Smooth 

70% 

Pref. 

EFM 

Smooth 

90% 

Pref. 

EFM 

Smooth 

95°/ 

Pref. 

>  EFM 
Smooth 

£,  [GPa] 

5.074 

4.642 

5.390 

4.914 

5.944 

5.453 

6.424 

5.923 

E 

■^23 

[GPa] 

5.836 

5.679 

6.217 

6.045 

6.818 

6.586 

7.282 

6.990 

G23  [GPa] 

0.700 

0.735 

0.773 

0.818 

0.839 

0.889 

0.870 

0.920 

G12 

[GPa] 

0.944 

0.957 

1.037 

1.050 

1.144 

1.157 

1.213 

1.226 

^21 

0.319 

0.357 

0.316 

0.355 

0.311 

0.350 

0.306 

0.344 

V32 

0.232 

0.226 

0.227 

0.221 

0.222 

0.219 

0.218 

0.218 
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The  preferential  fiber  orientation  exhibits  greater  stiffness  for  all  %  EFM  values 
tested.  Additionally,  the  stiffnesses  and  shear  moduli  increase  for  increasing  %  EFM,  as 
is  expected.  The  shear  moduli  in  the  23  direction  and  Poisson’s  ratio  in  the  23  direction 
are  similarly  matched  between  models. 

Table  23.  Cortical  bone  in  tension. 


50% 

Pref. 

EFM 

Smooth 

70% 

Pref. 

EFM 

Smooth 

90% 

Pref. 

EFM 

Smooth 

95% 

Pref. 

EFM 

Smooth 

El  [GPa] 

4.110 

3.704 

4.396 

3.930 

4.622 

4.120 

4.696 

4.193 

E 

^23 

[GPa] 

4.779 

4.649 

5.138 

5.009 

5.408 

5.286 

5.485 

5.371 

G,,[GPa] 

0.666 

0.699 

0.735 

0.778 

0.798 

0.846 

0.829 

0.876 

Gn 

[GPa] 

0.887 

0.900 

0.970 

0.983 

1.039 

1.054 

1.069 

1.085 

^21 

0.375 

0.427 

0.378 

0.432 

0.398 

0.455 

0.423 

0.482 

*32 

0.257 

0.250 

0.256 

0.247 

0.268 

0.257 

0.284 

0.272 

The  results  for  tension  are  similar  to  those  in  compression.  The  preferential  model 
exhibits  greater  stiffnesses  for  all  models.  Additionally,  the  shear  moduli  and  Poisson’s 
ratio  in  the  23  direction  are  also  similarly  matched  for  models.  For  both  compressive  and 
tensile  tests,  the  preferential  and  smooth  models  produce  similar  results.  For  this  reason, 
we  can  assume  that  the  physical  orientation  of  the  fibers  does  not  drastically  affect  the 
stiffness  of  cortical  bone. 

3.  Cortical  Model  Discussion 

The  results  of  the  cortical  model  can  be  compared  to  results  from  experimental 
and  theoretical  testing.  These  results  are  shown  in  Table  24. 


Table  24.  Comparison  of  values  for  cortical  results. 


Elastic  Modulus 

Method 

Source 

[GPal 

6 

FEM,  20%  Mineral  Volume  Fraction 

[1] 

56 


Elastic  Modulus 
[GPa] 

Method 

Source 

9 

FEM,  30%  Mineral  Volume  Fraction 

[1] 

15 

FEM,  40%  Mineral  Volume  Fraction 

[1] 

21 

FEM,  50%  Mineral  Volume  Fraction 

[1] 

7.28 

Microdisplacements  by  Machine 
Vision  Photogrammetry 

[73] 

17.5±1.9 

Uniaxial  Tension 

[74] 

17.8±2.1 

Uniaxial  Tension 

[75] 

12-22 

Assumptions 

[8] 

19.1±5.4 

Nanoindentation  (Transverse) 

[76] 

17.1±3.15 

Uniaxial  Compression  and  Tension1 

[77] 

17.0 

Uniaxial  Compression  and  Tension 

[78] 

14.91±0.52 

Acoustic  Microscopy  (Transverse) 

[79] 

20.55±0.21 

Acoustic  Microscopy  (Longitudinal) 

[79] 

16.58±0.32 

Nanoindentation  (Transverse) 

[79] 

23.45±0.21 

Nanoindentation  (Longitudinal) 

[79] 

22.5±1.3 

Nanoindentation  (Longitudinal) 

[80] 

1  Only  two  femurs  sampled  had  statistically  significant  differences  in  compression  and 
tension 


The  experimental  stiffnesses  shown  in  Table  24  determine  that  our  model 
underestimates  the  longitudinal  stiffness  of  cortical  bone.  This  could  be  due  to  several 
factors,  but  it  is  most  likely  due  to  the  assumptions  made  at  various  hierarchies. 
Furthermore,  the  similarity  of  tensile  and  compressive  testing  of  human  cortical  samples 
could  be  due  to  the  posthumous  treatment  of  samples.  Bone  is  kept  wet  throughout 
machining  and  processing,  however,  saline  solution  is  not  always  used.  The  use  of 
regular  water  can  affect  the  ion  balance  and  presence  of  water  at  the  nanoscale. 
Additionally,  the  bone  is  not  always  kept  at  physiological  temperatures.  The  most 
profound  effects  may  be  due  to  freezing  of  the  samples.  The  expansion  of  the  liquid  in 
the  void  spaces  may  have  altered  the  results.  Additionally,  water  exists  as  a  very  tightly 
bound  column  in  TC  molecules.  The  freezing  of  these  molecules  may  disrupt  their 
structure,  affecting  results  in  tension  and  compression. 

The  finite  element  model  listed  for  comparison  in  Table  24  found  lower  estimates 
of  stiffnesses  for  two  of  the  four  mineral  volume  fractions  tested.  Only  at  higher  mineral 
volumes  does  the  FEM  approach  the  values  found  through  experimental  testing.  The 


57 


stiffnesses  found  from  this  study  are  proportional  to  mineral  fraction  volume.  Table  19 
lists  the  mineral  volume  fractions  for  the  different  %EFM  tested.  This  study  utilizes  a 
maximum  of  25%  mineral  volume  fraction  which  yielded  a  longitudinal  stiffness  of 
approximately  7  GPa.  From  the  FEM  listed  above,  the  results  would  predict  a 
longitudinal  stiffness  of  approximately  7.5  GPa.  This  shows  that  the  relatively  simple 
models  presented  match  those  of  complex  finite  element  models.  It  can  also  be  reasoned 
that  increasing  the  mineral  volume  fraction  in  this  study  will  alter  the  macroscale  results. 

While  most  accepted  results  of  human  bone  material  properties  depict 
longitudinal  or  transverse  stiffness,  several  studies  have  been  completed  that  analyze  the 
orthotropic  properties  of  bone.  These  values  are  shown  in  Table  25. 


Table  25.  Elastic  moduli  of  compact  bone.  All  stiffnesses  are  in  GPa. 


Location 

Femur 

Tibia 

Femur 

Femur 

Method 

Ultrasound 

Mechanical 

Testing 

Ultrasound 

Ultrasound 

Preferential  Model 
95%EFM,  Compression 

Source 

[781 

[811 

[821 

[831 

This  Study 

En 

11.17 

6.94 

13.18 

18.57 

6.42 

E22 

11.17 

8.56 

14.56 

19.43 

6.39 

E* 

17.21 

18.45 

21.67 

28.29 

8.46 

G23 

3.3 

4.91 

6.56 

8.71 

0.87 

Gn 

3.3 

3.56 

5.85 

8.71 

0.97 

Gn 

3.6 

2.41 

4.74 

8.71 

1.62 

vn 

0.595 

0.495 

0.377 

0.323 

0.419 

V13 

0.298 

0.142 

0.237 

0.209 

0.191 

V23 

0.298 

0.119 

0.225 

0.126 

0.193 

The  values  depicted  above  are  the  orthotropic  material  properties.  The  raw  values 
of  the  preferential  fiber  model,  with  90%  EFM  and  in  compression,  are  shown  in  Table 
25.  These  values  do  not  compensate  for  transverse  isotropic  assumptions.  The  model 
presented  in  this  study  exhibits  similar  trends  as  the  experimental  tests.  The  stiffness  in 
the  3  direction  is  greater  than  both  the  stiffnesses  in  the  1  and  2  directions.  The 
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orthotropic  bone  properties  confirm  that  this  study  produces  an  underestimate  of  material 
properties.  Moreover,  the  results  shown  in  Table  25  oppose  the  assumption  that  small 
sections  of  macroscale  bone  are  transversely  isotropic. 

One  of  the  four  measurements  of  the  elastic  moduli  assumes  transverse  isotropy. 
The  other  three  measured  orthotropic  material  properties.  These  properties  are  defined  for 
small  samples,  not  single  osteons.  Additionally,  the  differences  between  transverse 
isotropic  and  orthotropic  measurements  are  small.  Large  deviations  have  been  credited  to 
regions  of  greater  porosity  [84],  As  a  whole,  cortical  bone  exhibits  orthotropic  symmetry, 
or  a  higher  degree  of  symmetry  [85].  This  means  that  the  minimum  level  of  symmetry  is 
orthotropic  symmetry,  but  regions  of  transverse  isotropy  or  isotropy  may  exist.  These 
higher  symmetries  can  be  attributed  to  the  packing  of  Haversian  systems  [85].  Hexagonal 
and  random  packing  can  be  identified  as  transversely  isotropic.  Due  to  the  preferential 
alignment  of  the  Haversian  canal  to  the  long  axis  of  the  bone,  the  assumption  that  cortical 
bone  is  transverse  isotropic  can  be  validated  for  ideal  materials.  The  imperfections  of 
living  tissues  may  cause  lower  order  symmetry. 

In  addition  to  the  anisotropic  properties  of  bone,  the  void  spaces  present  in  bone 
may  not  be  randomly  distributed.  The  ability  of  bone  to  sense  the  stress  field  is 
hypothesized  to  be  from  the  cell  processes  located  within  the  void  space  of  lamellar 
layers  [12].  It  can  be  extrapolated  that  the  canaliculi  are  aptly  located  so  as  to  optimally 
sense  the  mechanical  state  of  the  surrounding  bone.  Thus,  it  can  be  derived  that  the  void 
space  of  bone  is  a  non-random  occurrence. 

The  alternating  orientation  of  lamellar  layers  had  little  impact  on  the  stiffness  of 
cortical  bone.  However,  this  twisting  is  thought  to  aid  in  resisting  crack  propagation  [64], 
This  can  help  increase  fracture  toughness  and  prevent  damage  from  propagating  from  one 
lamellar  layer  to  the  next. 

C.  CANCELLOUS  BONE 

Cancellous  bone  is  a  porous  macrostructure  that  consists  of  small  plates  and  rods. 
The  porous  nature  of  cancellous  bone  allows  for  it  to  house  bone  marrow  and  to  reduce 
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the  weight  of  the  skeletal  system  [6],  The  individual  rods  of  cancellous  bone  are  known 
as  trabeculae. 

1.  Three  Dimensional  Cancellous  Structure 

Cancellous  bone  is  more  metabolically  active  than  cortical  bone.  This  activity 
leads  to  more  frequent  absorption  and  deposition  of  bone.  Additionally,  the  deposition  of 
new  bone  is  related  to  the  applied  stresses  on  the  cancellous  bone  [9],  [86],  [87].  As  a 
result  of  the  directed  deposition  of  lamellar  layers,  the  outermost  layers  of  trabecular 
bone  are  oriented  parallel  to  the  longitudinal  axes  of  the  trabeculae  [9].  The  adjacent 
layers  vary  slightly  in  their  orientation,  as  they  were  formed  under  a  different  state  of 
stress.  The  resulting  structure  of  a  trabecula  is  that  of  lamellar  layers  that  follow  a 
relatively  uniform  fiber  direction.  This  is  shown  in  Figure  19. 


Figure  19.  Trabecular  bone  lamellar  layer  orientation,  from  [9]. 

The  layers  intersect  at  low  angles,  with  respect  to  the  longitudinal  axis  of  the 


trabecula.  As  a  result,  the  properties  of  a  single  trabecula  can  be  found  through  a 
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laminated  composite  constructed  of  a  -10°/-57075°/10°  layered  motif.  The  center  of  the 
strut,  indicated  by  the  star  in  Figure  19,  shows  a  remnant  of  the  original  trabecular 
network  [9],  which  has  yet  to  be  regrown.  These  remnants  will  be  assumed  to  exhibit 
similar  material  properties  to  the  lamellar  layers.  This  will  cause  the  model  to  have  a 
slightly  overestimation  of  the  longitudinal  modulus,  since  more  lamellar  layers  are 
assigned  to  the  longitudinal  axis. 

While  a  laminated  composite  material  represents  the  material  properties  of  a 
single  trabecula,  it  does  not  quantify  the  properties  of  cancellous  bone.  The  structure  of 
cancellous  bone  is  a  uniquely  three  dimensional  problem,  that  cannot  be  solved  through  a 
two  dimensional  approximation  [88],  Additionally,  there  is  much  heterogeneity  in 
cancellous  bone  at  different  anatomical  locations  [89]. 

Early  models  of  cancellous  bone  utilized  a  model  of  rods  and  plates.  Four  early 
models  were  an  asymmetric  rod-like  cubic  model,  a  plate-like  cubic  model,  a  rod-like 
hexagonal  columnar  model,  and  a  plate-like  hexagonal  columnar  model  [87].  These  early 
models  helped  to  shape  simple  models  of  cancellous  bone.  However,  these  models 
provided  asymmetrical  properties  for  cancellous  bone.  Additionally,  the  early  models 
included  Euler  buckling  as  a  failure  mechanism  [87].  Later  experiments  have  shown  that 
failure  of  cancellous  bone  is  most  commonly  due  to  microscopic  cracking,  which 
removes  buckling  as  a  failure  mode  for  trabeculae  [90]. 

With  the  decrease  in  cost  and  increase  in  speed  of  computing,  complex  finite 
element  models  have  been  used  to  model  cancellous  bone.  Three  dimensional  models 
formed  from  micro-computed  tomography  can  replicate  small  sections  of  bone  [91]. 
Furthermore,  two  unit  cells  have  been  proposed  that  are  able  to  accurately  model 
cancellous  bone.  Kadir  et  al.  [91]  compared  the  results  of  prismatic  unit  cells  and 
tetrakaidecahedral  unit  cells  to  those  of  a  micro-computed  tomography  model.  The 
authors  found  that  both  unit  cells  accurately  represent  the  mechanical  properties. 
Additionally,  Guo  and  Kim  [92]  have  shown  that  a  complex  finite  element  model  of 
several  tetrakaidecahedron  cells  can  accurately  represent  different  levels  of  bone  loss  due 
to  ageing. 
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This  study  will  make  use  of  a  tetrakaidecakedron  unit  cell  to  calculate  the 
macroscale  properties  of  cancellous  bone.  Several  values  of  trabecular  bone  properties 
will  be  used,  as  well  as  different  bone  densities.  The  results  of  the  macroscale  properties 
will  then  be  compared  to  experimental  data. 

2.  Trabecular  Bone  Model 

Trabecular  bone  will  be  modeled  by  the  same  layered  composite  model  as  utilized 
for  cortical  bone.  However,  the  fiber  orientations  and  volume  percent  of  each  layer  will 
be  altered. 


a.  Trabecular  Model  Par  am  eters 

The  parameters  for  trabecular  bone  are  listed  in  Table  26. 


Table  26.  Trabecular  bone  layered  composite  model  parameters. 


Fiber 

Direction 

Volume 

Percent 

-10° 

20 

-5° 

20 

0° 

20 

5° 

20 

10° 

20 

Due  to  the  lack  of  Haversian  canals  in  trabecular  bone,  the  void  space  is  direct  a 
result  of  microscopic  canaliculi.  This  void  space  will  be  evaluated  at  5%  volume  percent 
and  will  be  treated  as  a  liquid  in  compression  and  a  void  in  tension. 

b.  Trabecular  Results 

The  results  of  the  trabecular  laminated  composite  model  for  both  compression  and 
tension  are  shown  in  Table  27. 
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Table  27.  Trabecular  bone  in  compression  and  tension. 


50% : 

Comp. 

EEM 

Tens. 

70% 

Comp. 

EFM 

Tens. 

90% : 
Comp. 

EFM 

Tens. 

95%  ] 
Comp. 

iFM 

Tens. 

Ex  [GPa] 

6.047 

5.099 

6.427 

5.158 

7.848 

5.223 

9.451 

5.245 

E2i  [GPa] 

8.519 

7.762 

9.186 

8.569 

9.782 

9.134 

9.975 

9.252 

^23  [GPa] 

1.259 

1.197 

1.461 

1.389 

1.605 

1.526 

1.636 

1.555 

Gn  [GPa] 

0.658 

0.623 

0.672 

0.636 

0.717 

0.673 

0.773 

0.720 

V2X 

0.242 

0.254 

0.235 

0.264 

0.236 

0.318 

0.241 

0.382 

V32 

0.234 

0.253 

0.232 

0.241 

0.240 

0.237 

0.248 

0.237 

The  results  in  compression  are  stronger  than  the  results  in  tension  for  all  models. 
The  difference  between  tensile  and  compressive  stiffnesses  magnifies  for  increasing 
%EFM.  The  results  of  the  trabecular  bone  model  can  be  compared  to  experimental  and 
theoretical  results.  Early  models  tested  macroscale  cancellous  bone  in  order  to  back 
calculate  the  material  properties  of  trabecular'  bone.  Current  methods  directly  measure 
trabecular  material  properties  through  the  use  of  microelectromechanical  systems  and 
nanoindentation.  Table  28  lists  the  experimental  and  theoretical  results  for  trabecular 
bone. 


Table  28.  Comparison  of  values  for  trabecular  results. 


Elastic  Modulus 
[GPa] 

Method 

Source 

18.5±1.5 

Assumption 

[93] 

11.38 

Inelastic  Buckling 

[94] 

12.7±2.0 

Ultrasound  (Isotropic) 

[95] 

15 

Microharduess  (Transverse) 

[96] 

3.81 

Three-point  bending 

[97] 

5.35±1.36 

Four -point  bending 

[98] 

10.4±3.5 

Tensile  Test 

[99] 
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14.8±1.4 

Ultrasoiuid  (Isotropic) 

[99] 

13.4±2.0 

Nanoindentation  (Transverse) 

[80] 

18.0±2.8 

Tension  and  Compression  Tests 

[75] 

1 1.4±5.6 

Nanoindentation  (Transverse) 

[76] 

17.5±1.12 

Acoustic  Microscopy  (Transverse) 

[79] 

18.14±1.7 

Nanoindentation  (Transverse) 

[79] 

The  results  for  trabecular  bone  are  similar  to  those  of  cortical  bone,  as  the  model 
used  in  this  study  underestimates  the  stiffness  values.  However,  the  trabecular  model  is 
closer  to  the  mean  values  of  accepted  longitudinal  and  transverse  results.  The  values 
established  from  bending  tests  were  found  by  testing  machined  samples  of  bone.  The 
same  post-processing  issues  that  affected  the  measurement  of  cortical  bone  will  affect  the 
measurements  of  trabecular  bone. 


3.  Tetrakaidecahedron  Finite  Element  Model 

The  computation  of  a  unit-cell  model  for  cellular  materials  was  developed  by 
Kwon,  et.  al.  [100],  This  model  utilizes  a  tetrakaidecahedron  (4KDH)  unit  cell  as  the  base 
unit.  The  4KDH  is  a  regular  truncated  octahedron  and  was  developed  by  Lord  Kelvin  in 
1887  as  an  attempt  to  model  soap  bubble  formation  in  foam.  A  4KDH  provides  minimum 
surface  area  for  a  given  volume  [100],  Additionally,  the  4KDH  cell  is  space  filling, 
meaning  that  its  faces  uniformly  intersect  with  adjoining  cells  [101],  Visualization  of  the 
4KDH  structure  is  shown  in  Figure  20. 


Figure  20.  Tetrakaidecahedron,  from  [101], 
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The  4KDH  has  eight  hexagonal  faces,  six  square  faces,  36  edges  and  24  vertices. 
For  the  open  cell  structures  present  in  cancellous  bone,  the  faces  are  treated  as  voids  and 
the  edges  are  treated  as  trabeculae.  A  simple  finite  element  model  was  developed  by 
[100]  to  analyze  the  4KDF1  frame. 

The  4KDF1  unit  cell  is  a  cube,  where  the  frame  structure  is  composed  of  uniform 
circular  beams.  The  cross-section  of  the  frame  elements  are  dependent  on  the  volume 
fraction  defined  for  the  unit-cell  [100].  Each  frame  element  is  defined  by  one  node  at 
each  end.  The  model  utilizes  six  degrees  of  freedom:  displacement  and  rotation  about 
each  of  the  three  axes.  Loads  applied  to  the  model  are  applied  at  the  vertices.  It  is 
assumed  that  the  beams  do  not  buckle.  Furthermore,  the  finite  element  model  developed 
by  [101]  was  developed  through  Galerkin’s  methods  and  linear  shape  functions,  and 
assumes  linear  elastic  behavior. 

The  finite  element  model  is  based  on  forces  applied  to  the  model  and  the  stiffness 
of  the  edge  elements.  This  relationship  is  defined  in  Equation  30. 

[A]{r/}={F}  (30) 


where  [A]  is  the  global  stiffness  matrix,  {<i}  is  the  nodal  displacement  and  {F}  is  the 
applied  forces.  The  elemental  stiffness  matrix  is  defined  in  Equation  31  [101]. 


K 

Kn 

K\x 

Ke21 

(31) 


Furthermore,  the  stiffness  of  each  frame  element  must  be  rotated  based  on  its  direction 
cosines  in  order  to  include  it  in  the  global  stiffness  matrix.  The  complete  matrix 
representation  of  the  4KDH  model  is  shown  in  the  Appendix. 

The  FEM  model  for  the  4KDH  produces  isotropic  material  properties.  Therefore, 
from  an  applied  load,  the  material’s  Young’s  modulus  can  be  calculated  from  the 
resultant  strains.  Additionally,  since  each  frame  element  is  treated  as  an  isotropic 
material,  two  variations  will  be  assessed:  one  using  the  compressive  material  properties 
and  one  using  the  tensile  material  properties.  This  will  create  a  bound  on  the  material 
properties  of  cancellous  bone  for  each  percent  extrafibrillar  mineralization  used. 
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4. 


Cancellous  Bone  Model  Parameters 


In  order  to  assess  the  material  properties  of  cancellous  bone,  several  parameters 
must  be  identified.  The  4KDH  model  calls  for  isotropic  material  properties.  Since  the 
frame  elements  of  the  model  are  long  beams,  the  material  properties  of  the  trabecular 
long  axis  will  be  the  dominant  material  properties.  Therefore,  the  longitudinal  stiffness 
and  shear  modulus  will  be  used.  These  properties  for  the  different  %EFM,  in  tension  and 
in  compression,  are  listed  in  Table  27.  Furthermore,  the  effects  of  bone  marrow  will  be 
ignored  in  this  study.  Bone  marrow  has  little  ability  to  resist  shear  [102],  diminishing  its 
ability  to  act  as  an  elastic  foundation  for  the  beam  elements. 

The  radius  of  the  beam  elements  and  the  volume  density  of  the  bone  must  be 
defined.  A  universally  accepted  set  of  indices  have  been  identified  that  help  to  address 
cancellous  structure:  bone  volume  per  tissue  volume  (BV/TV),  bone  sirrface  per  tissue 
volume  (BS/TV),  bone  sirrface  per  bone  volume  (BS/BV)  [89],  These  indices  can  help  to 
derive  the  structure  of  the  trabeculae:  trabecular  thickness  (Tb.Th),  trabecrrlar  separation 
(Tb.Sp),  and  trabecular  number  (Tb.N)  [89]. 

a.  Model  Overview 

Cancellous  indices  for  various  anatomical  locations  have  been  defined.  Table  29 
lists  the  two  indices  most  pertinent  to  the  4KDH  model  and  identifies  the  anatomical 
location  of  the  measurements. 


Table  29.  Anatomical  cancellous  bone  indices. 


BV/TV 

Tb.Th  (mm) 

Anatomical  Location 

Source 

0.2654 

Femoral  Neck 

[103] 

0.2067 

0.172 

Femoral  Head 

[104] 

0.27 

0.2 

Distal  Femur 

[105] 

0.27 

0.19 

Proximal  Femur 

[105] 

0.105 

Greater  Trochanter 

[103] 

0.111 

Tibia 

[103] 

0.34 

0.36 

Radius 

[106] 

0.41 

0.33 

Distal  Tibia 

[107] 

0.27 

0.47 

Distal  Tibia 

[107] 

0.186 

Average  (Young  Age)1 

[108] 

0.174 

Average  (Medium  Age)1 

[108] 
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Young  age  is  defined  as  16-39  years-old  and  medium  age  is  defined  as  40-59  yeais-old 


The  measurements  taken  by  [107]  were  used  for  comparison  of  high-resolution 
peripheral  quantitative  X-ray  microtomography  to  the  standardized  measurements  taken 
by  magnetic  resonance  imaging.  Therefore,  the  differences  in  BV/TV  and  Tb.Th  of  the 
distal  tibia  are  due  to  the  different  measurement  techniques.  This  study  will  make  use  of 
the  measurements  of  the  femur,  as  it  is  the  most  tested  bone  in  the  human  anatomy. 

b.  Tetrakaidecahedron  Model  Parameters 

In  order  to  utilize  the  4KDH  model,  the  diameters  of  the  beam  elements  are  set  to 
the  average  thickness  of  the  femoral  trabeculae.  This  was  calculated  as  0.19  mm,  from 
the  indices  defined  in  the  femoral  measurements.  Additionally,  the  BV/TV  index  was 
used  to  define  the  beam  element  length.  The  length  was  calculated  so  that  the  relative 
volume  of  the  4KDH  unit  cell  matches  that  of  the  trabecular  unit  cell.  The  load  applied  to 
each  model  was  derived  from  the  area  of  the  4KDH  unit  cell  face  so  that  a  constant  stress 
was  analyzed.  This  study  will  calculate  the  effective  stiffnesses  of  the  three  femoral 
indices:  the  femoral  neck,  the  femoral  head  and  the  distal/proximal  femur.  All 
calculations  assume  a  constant  trabecular  diameter  of0.19  mm. 

c.  Model  Results 

The  results  of  the  4KDH  model  were  calculated  using  the  trabecular  results  in 
compression  and  tension.  These  two  values  will  be  used  to  create  bounds  on  the  stiffness 
of  cancellous  bone.  The  results  of  the  three  anatomical  locations  are  shown  in  Table  30. 


Table  30.  Cancellous  model  results,  all  stiffnesses  expressed  in  MPa. 


50%  I 
Comp. 

EFM 

Tens. 

70%  I 
Comp. 

ZFM 

Tens. 

90%  1 
Comp. 

ZFM 

Tens. 

95%  1 
Comp. 

EFM 

Tens. 

Femoral  Neck 

292 

246 

311 

249 

379 

252 

457 

253 

Femoral  Head 

196 

165 

209 

167 

255 

170 

307 

170 

Distal/Proximal 

300 

253 

319 

256 

390 

259 

469 

260 
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The  results  of  the  cancellous  model  are  stiffer  in  compression  than  tension  for  all 
models.  Additionally,  the  distal  and  proximal  femur  is  stiffer  than  both  the  femoral  neck 
and  head.  This  is  expected,  due  to  its  greater  density.  The  results  of  the  cancellous  model 
can  be  compared  to  experimental  and  theoretical  results. 

5.  Cancellous  Model  Discussion 

The  accepted  values  of  cancellous  bone  are  shown  in  Table  3 1 . 


Table  3 1 .  Comparison  of  values  for  cancellous  results. 


Elastic  Modulus 
[MPa] 

Method 

Source 

47±2 

Finite  Element  Method 

[109] 

300±200 

Micro-compression  testing 

[110] 

1 10±90 

Compressive  Testing 

[111] 

385±315 

Estimation 

[8] 

550±450 

Mechanical  Testing 

[112] 

The  accepted  values  are  in  agreement  with  the  values  found  in  this  study. 
However,  due  to  the  underestimation  of  previous  hierarchies,  a  similar  underestimate 
should  be  expected  for  cancellous  bone.  The  results  of  the  cancellous  model  produce 
elastic  moduli  greater  than  expected.  This  could  be  due  to  several  causes.  The  assitmed 
isotropy  of  cancellous  bone  is  a  broad  assumption.  Orthotropic  elastic  moduli  have  been 
measured  through  ultrasonic  measurements  [113].  However,  the  calculation  of 
orthotropic  materials  tor  this  study  would  require  a  complex  model  based  on  a  specific 
anatomical  site.  Utilizing  a  generalized  4KDH  unit  cell,  an  isotropic  stiffness  was 
calculated  that  confirms  the  validity  of  previous  hierarchies.  This  is  a  computationally 
inexpensive  model  that  is  not  site  specific. 

The  results  of  the  cancellous  model  help  identify  the  reduced  stiffness  of  the 
porous  structure.  The  physiological  purpose  of  cancellous  bone  is  to  provide  lightweight 
structural  support  to  abnormal  volumes  of  skeletal  bone.  Additionally,  the  porous 
structure  of  cancellous  bone  provides  a  cellular  infrastructure  for  bone  marrow,  creating  a 
large  surface  area  for  metabolic  activity.  Furthermore,  the  in  vivo  trabecular  network  has 
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been  optimized  through  repeated  applied  stresses.  The  trabeculae  respond  to  the  applied 
mechanical  loads  [9],  [86],  [87].  Due  to  this,  the  in  vivo  mechanical  properties  of 
cancellous  bone  should  be  greater  than  those  found  through  experimental  results. 
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V.  ANALYSIS 


In  order  to  analyze  the  robust  nature  of  the  models  presented  in  this  study, 
changes  to  the  different  hierarchies  will  be  assessed.  Additionally,  the  key  hierarchy  will 
be  identified  that  has  the  greatest  influence  on  the  macroscale  values.  In  order  to  assess 
the  overall  effects  of  the  model,  we  must  establish  a  control  model.  From  the  results  of 
the  fiber  model,  we  can  expect  a  high  level  of  extrafibrillar  mineralization  within  living 
bone.  Additionally,  the  linear  fibril  model  was  assumed  to  be  the  proper  fibril  crystal 
packing  method  due  to  the  assessed  mineral  volume.  Therefore,  the  control  model  will  be 
comprised  of  the  linear  fibril  HA  packing  and  95%  EFM.  Moreover,  the  model  will  only 
be  tested  in  compression.  The  macroscale  cortical  model  will  assume  a  preferential  fiber 
orientation,  since  the  results  of  the  two  cortical  models  are  similar. 

The  results  that  will  be  assessed  for  macroscale  models  will  be  longitudinal 
stiffness,  transverse  stiffness,  and  mineral  volume  fraction.  Additionally,  the  results  of 
the  modified  hierarchy  will  be  assessed  in  comparison  to  cortical  and  trabecular  values, 
as  the  results  of  cancellous  bone  are  not  as  robust  as  those  of  trabecular  bone.  A  separate 
analysis  will  be  completed  to  assess  bone  loss.  In  order  to  analyze  the  effects  of  bone  loss 
on  the  properties  of  cancellous  bone,  different  bone  densities  will  be  evaluated  using  the 
FEM  presented  for  macroscale  cancellous  bone. 

A.  MODIFIED  HIERARCHY 

This  analysis  will  make  changes  at  all  microscale  and  macroscale  hierarchies  in 
order  to  assess  their  contribution  to  cortical  and  trabecular  properties.  At  the  microscale 
level,  the  fibrillar  subunit  will  be  altered.  Additionally,  the  diameter  of  the  fibrils  and  the 
thickness  of  the  EFM  will  be  adjusted.  At  the  macroscale,  the  dimensions  of  the  fiber  and 
the  disordered  matrix  will  be  adjusted.  Lastly,  both  the  Voigt  and  Reuss  averages  of 
cortical  and  trabecular  bone  will  be  found.  These  two  models  create  an  upper  and  lower 
bound  on  the  stiffness  of  a  composite  based  on  volume  fraction  of  each  constituent  [114]. 
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1.  Hierarchical  Adjustments 

The  adjustments  to  the  hierarchy  will  be  completed  individually.  The  adjustments 
will  fir  st  be  described  and  the  results  of  all  adjustments  will  be  compared  to  the  control 
model.  The  first  alteration  will  be  made  to  the  unit  cell  of  the  fibrillar  model.  The 
dimensions  of  the  HA  ciystal  will  remain  the  same,  but  the  volume  of  the  surrounding  TC 
matrix  will  be  reduced.  The  dimensions  shown  in  Table  32  will  be  used. 

Table  32.  Unit  cell  alteration  to  fibrillar  model. 


Dimension 

Model  1 
Value  [run] 

Model  2 
Value  [ mill 

Model  3 
Value  [nml 

a\ 

50 

50 

50 

a2 

17 

17 

17 

25 

25 

25 

3 

1 

1 

cx 

3 

3 

3 

ci 

6 

6 

3 

These  alterations  will  increase  the  mineral  volume  percent  of  the  fibrillar  model 
while  also  changing  the  material  properties.  The  next  hierarchy  to  be  altered  is  that  of  the 
bone  fiber.  The  model  of  the  fiber  relies  on  volume  percent  of  both  the  bone  fibril  and  the 
extra  fibrillar  mineral.  Two  variations  will  be  tested  that  alter  the  volume  percent  of  the 
fibril.  The  original  model  assumed  a  72.25  vol%  fibril.  The  first  alteration  will  assume  a 
62.25  vol%  fibril  and  the  second  alteration  will  assume  an  82.25  vol%  fibril.  These 
variations  will  be  referred  to  as  Model  4  and  Model  5,  respectively. 

The  only  component  of  the  macroscale  hierarchies  that  will  be  adjusted  is  that  of 
lamellar  bone.  The  properties  of  lamellar  bone  were  calculated  by  discerning  the  volume 
percent  of  bone  fibers  in  the  disordered  fibrillar  matrix.  The  original  model  was 
composed  of  only  38.44  vol%  fiber.  The  adjusted  models  will  increase  the  fiber  volume 
percent  in  order  to  compare  results.  Model  6  will  assume  50  vol%  fiber.  Model  7  will 
assume  70  vol%  fiber,  and  Model  8  will  assume  90  vol%  fiber. 
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2.  Adjustment  Results 

The  results  of  Models  1  through  8  are  shown  in  Table  33.  These  models  adjust 
only  the  hierarchy  discussed,  leaving  all  other  levels  with  then  original  geometry. 


Table  33.  Results  of  adjusted  hierarchies,  stiffnesses  shown  in  GPa. 


Model 

Vo 

Cortica 

igt 

E 

J~'23 

/  Bone 
Re 

*r 

uss 

*23 

Vc 

Trabecu 

igt 

E 

-^23 

lar  Bom 
Re 

*r 

0 

uss 

E 

-^23 

Mineral 
Volume 
Fraction  [%] 

1 

7.05 

8.03 

3.75 

4.48 

10.28 

11.19 

8.04 

9.00 

29.97 

2 

8.27 

9.53 

3.92 

4.72 

11.63 

13.87 

8.76 

10.48 

31.50 

3 

10.26 

11.87 

4.33 

5.15 

14.18 

17.73 

10.20 

12.49 

42.18 

4 

6.95 

7.87 

3.58 

4.33 

10.79 

10.64 

8.15 

8.60 

28.03 

5 

5.84 

6.61 

3.61 

4.23 

8.15 

9.07 

6.84 

7.72 

22.00 

6 

7.21 

8.30 

3.35 

4.21 

10.49 

11.84 

7.77 

9.12 

27.53 

7 

8.93 

10.49 

2.75 

3.81 

12.29 

16.07 

7.68 

10.34 

31.89 

8 

11.25 

13.39 

1.81 

2.84 

14.08 

22.04 

6.38 

10.05 

36.22 

Original 

6.42 

7.28 

3.58 

4.29 

9.45 

9.98 

7.52 

8.24 

25.01 

The  results  show  that  the  transverse  stiffness  is  greater  than  longitudinal  stiffness 
for  all  models.  A  general  trend  can  be  seen  that  increasing  mineral  volume  fraction 
increases  macroscale  stiffnesses.  Additionally,  reducing  the  TC  matrix  within  the  fibrillar 
subunit  increases  macroscale  properties.  Increasing  the  volume  percent  of  the  fibril 
within  the  fiber  model  has  an  inverse  effect  on  its  stiffness.  This  is  due  to  the  extra 
fibrillar  mineral.  The  less  volume  that  is  occupied  by  the  fibril,  the  more  volume  is 
occupied  by  the  mineral.  Increasing  the  size  of  the  bone  fibers  within  the  lamellar  layers 
increases  the  macroscale  stiffnesses. 

Furthermore,  several  trends  can  be  identified  between  the  Voigt  and  Reuss 
models.  The  Voigt  model  assumes  isostrain  conditions,  whereas  the  Reuss  model 
assumes  isostress  conditions  [114].  Thus,  the  Voigt  is  the  upper  bound  of  the  material 
properties  while  the  Reuss  is  the  lower  bound.  These  are  derived  by  averaging  the 
stiffness  matrices  or  the  compliance  matrices,  respectively.  The  results  of  the  cortical 
model  show  a  large  difference  between  the  Voigt  and  Reuss  models.  This  is  due  to  the 
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large  theta  values  of  the  rotated  layers.  The  trabecular  models  show  less  variation  because 
the  layers  are  constrained  to  small  thetas. 

Furthermore,  the  most  influential  hierarchy  is  the  fibrillar  model.  Small  changes 
in  the  basic  hierarchy  are  compounded  throughout  the  upper  levels.  Changes  to  the 
fibrillar  model  can  be  physically  embodied  by  reducing  the  vertical  spacing  between  HA 
crystals.  The  second  most  influential  hierarchy  is  that  of  the  lamellar  layers.  The 
relatively  small  volume  percent  of  fibers  pr  esent  in  the  original  model  allowed  for  large 
increases  in  fiber  diameter.  This  increase  in  fiber  diameter  decreases  the  relative 
proportion  of  fibrillar  matrix. 

3.  Optimized  Adjustment 

The  parameters  of  several  hierarchies  were  adjusted  simultaneously  in  an  attempt 
to  optimize  the  macroscale  results.  At  the  microscale,  both  the  fibrillar  subunit  and  fiber 
models  were  adjusted.  The  fibrillar  subunit  utilized  the  unit  cell  dimensions  of  Model  2. 
The  fiber  model  assumed  a  fibril  vol%  of  60%.  The  lamellar  model  assumed  a  fiber  vol% 
of  70%.  The  results  from  the  model  are  shown  in  Table  34. 


Table  34.  Optimized  macroscale  bone  properties. 


Coi 

Voigt 

tical 

Reuss 

Trab 

Voigt 

ecular 

Reuss 

Ex  [GPa] 

12.71 

2.82 

16.72 

8.91 

E„  [GPa] 

15.03 

3.98 

24.17 

12.73 

G23  [GPa] 

0.914 

0.004 

2.18 

0.020 

Gn  [GPa] 

1.62 

0.004 

0.52 

0.019 

V2l 

0.272 

0.413 

0.165 

0.304 

V32 

0.153 

0.176 

0.164 

0.176 

The  optimized  hierarchy  of  bone  contained  43.82%  mineral  volume  fraction.  The 
values  of  Ex  for  both  models  and  the  value  of  E2 3  for  the  cortical  model  are  within  the 
range  of  accepted  values.  Transverse  stiffness  of  the  trabecular  model  is  gr  eater  than 


74 


expected.  Additionally,  the  results  can  be  compared  to  the  orthotropic  material  properties 
of  cortical  bone  shown  in  Table  25.  The  cortical  shear  moduli  and  vn  values  are  slightly 
lower  than  the  tested  values,  although  the  Young’s  moduli  are  within  normal  bounds. 

The  accepted  stiffnesses  of  cortical  bone  are  greater  than  those  of  trabecular  bone. 
The  results  from  this  study  show  that  the  trabecular  material  is  stiffer  than  the  cortical 
material.  Current  theories  accept  that  the  only  difference  between  trabecular  and  cortical 
material  is  the  volume  fraction  of  fiber  orientations,  which  is  respected  in  this  study. 
Therefore,  discrepancies  in  material  stiffness  are  due  to  variations  at  the  microscale 
hierarchies. 

This  could  be  due  to  several  assumptions,  such  as  the  assumption  that  the  lamellar 
layers  at  the  center  of  trabecular  stmts  exhibit  the  same  material  properties  as  the  younger 
bone.  The  layers  should  physically  portray  lamellar  layers  that  are  similar  in  orientation 
to  the  surface  layers.  However,  this  assumption  is  unique  for  each  individual  stmt  based 
on  its  historical  loading.  It  can  be  hypothesized  that  each  trabecular  stmt  exhibits 
different  mechanical  properties,  independent  of  its  dimensions.  This  makes  the 
quantification  of  cancellous  material  properties  very  difficult,  even  with  complex  finite 
element  models. 

Additionally,  the  differences  could  be  due  to  the  simplification  of  void  space 
within  the  macrostmctures.  The  cortical  bone  models  account  for  25%  void  space  due  to 
the  combination  of  Haversian  canals  and  canaliculi.  The  trabecular  model  identifies  only 
5%  void  space  due  to  the  lack  of  Haversian  canals.  However,  the  void  spaces  were 
theorized  to  be  in  locations  that  improved  the  mechanosensation  of  the  cells  that  directed 
bone  growth.  Therefore,  the  complex  geometry  of  cortical  and  trabecular  bone  may 
require  more  than  a  volume  fraction  analysis  to  determine  the  implications  of  void  space. 

B.  EFFECTS  OF  BONE  LOSS 

Bone  loss  is  an  identifiable  side  effect  of  both  age  and  disease.  Several  studies 
have  identified  quantitative  differences  between  healthy  individuals  and  individuals 
suffering  from  bone  loss.  The  differences  are  shown  with  respect  to  the  cancellous 

indices  of  BV/TV  and  Tb.Th.  The  effects  of  bone  loss  are  shown  in  Table  35. 
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Table  35.  Cancellous  bone  indices  for  individuals  with  bone  loss. 


BV/TV 

Tb.Th 

(rum) 

Anatomical 

Location 

Group  Identifier 

Source 

0.26 

0.33 

Distal  Radius 

Osteoporotic  Postmenopausal  Females 

[106] 

0.157 

Proximal  Tibia 

60-79  years-old 

[108] 

0.16 

Iliac  Crest 

70-90  years-old 

[115] 

While  no  results  are  shown  for  the  human  femur,  the  indices  shown  in  Table  35 
identify  both  a  reduction  in  BV/TV  and  a  thinning  of  trabeculae.  In  order  to  analyze  the 
effects  of  reductions  in  bone  density,  both  modes  of  bone  loss  will  be  tested.  The 
cancellous  model  will  assume  baseline  values  of  0.27  and  0.19  mm  for  BV/TV  and 
Tb.Th,  respectively.  The  models  will  utilize  the  optimized  material  properties  of 
trabecular  bone,  as  found  through  Voigt  averaging,  and  will  test  only  in  compression.  The 
baseline  cancellous  results  are  re-calculated  for  the  optimized  trabecular  properties. 


Table  36.  Cancellous  bone  material  properties  after  bone  loss. 


BV/TV 

Tb.Th 

[mm] 

Stiffness 

[MPa] 

0.27 

0.19 

830 

0.24 

0.19 

689 

0.20 

0.19 

514 

0.16 

0.19 

357 

0.27 

0.17 

830 

0.27 

0.15 

830 

0.27 

0.13 

830 

The  results  of  the  cancellous  model  highlight  a  weakness  of  the  finite  element 
model.  The  FEM  utilizes  a  ligament  radius  and  length.  In  order  to  assign  a  BV/TV  ratio, 
the  radius  is  defined  and  the  ligament  length  is  calculated.  Changes  to  the  trabecular 
thickness  do  not  affect  cancellous  stiffness  as  the  length  of  the  elements  is  also  reduced  to 
account  for  a  steady  BV/TV. 

The  results  show  a  steady  decrease  in  macroscale  stiffness  as  a  result  of 
decreasing  density.  This  shows  that  the  effects  of  a  10%  reduction  in  density  results  in  a 
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50%  reduction  in  stiffness.  Furthermore,  the  effects  of  the  trabecular  thinning  may  apply 
more  to  the  failure  modes  of  trabecular  bone  than  to  the  material  properties. 


C.  HELICAL  MODEL 

Early  models  of  osteons  speculated  that  the  angular  orientation  of  bone  fibers 
acted  as  spring  stiffeners.  They  assumed  that  the  twisting  fibers  within  a  lamellar  layer 
act  as  springs  [9].  The  stiffness  of  an  osteon  can  then  be  found  by  summing  the 
stiffnesses  of  the  individual  springs.  This  is  similar  to  the  model  used  for  the  calculation 
of  tropocollagen  material  properties.  We  can  test  the  validity  of  this  assumption  by 
calculating  the  spring  constant  of  a  single  bone  osteon.  This  can  be  done  by  summing  the 
spring  constants  of  all  the  fibers  within  an  osteon  using  Equation  1,  which  is  reproduced 
below. 


F  _  Gd4 

y  ~  8  d3na 


(i) 


The  shear  modulus  from  the  optimized  hierarchy  is  used.  The  diameter  of  the 
spring  will  vary  based  on  the  lamellar  layer.  This  model  will  assume  a  lamellar  thickness 
of  2.5  jum.  The  number  of  coils  will  depend  on  the  fiber  orientation  and  the  height  of  an 
osteon.  We  will  assume  half  of  the  fibers  are  oriented  at  30°,  half  of  the  fibers  are 
oriented  at  60°,  and  that  an  osteon  is  5  mm  tall.  This  fiber  orientation  is  similar  to  that  of 
the  preferential  fiber  model  presented  for  cortical  bone.  Additionally,  the  number  of 
fibers  contained  within  each  lamella  is  defined  by  the  circumference  of  the  layer  and  the 
diameter  of  the  fiber.  We  will  assume  a  fiber  diameter  of  2.10  /um ,  which  was  derived 
from  the  optimized  model.  We  can  then  find  the  spring  constant  of  each  lamellar  layer  by 
multiplying  the  individual  spring  constant  of  a  single  fiber  by  the  number  of  fibers  within 
the  layer.  The  total  spring  constant  of  the  osteon  can  be  found  by  summing  the  spring 
constants  of  the  lamellar  layers.  The  Young’s  modulus  of  the  osteon  can  then  be  found  by 
applying  Equation  2.  The  resulting  longitudinal  Young’s  modulus  is  found  to  be 
approximately  170  KPa. 
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This  model  produces  very  small  values  for  the  expected  stiffness  of  cortical  bone. 
The  results  are  more  than  two  orders  of  magnitude  less  that  what  we  would  expect.  For 
this  reason,  these  early  models  are  discredited. 
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VI.  ADVANCED  TOPICS 


The  models  presented  in  this  study  analyze  the  hierarchical  relationship  of  bone. 
This  chapter  will  serve  as  an  expansion  of  the  models  presented  earlier.  Additionally,  this 
section  discusses  several  topics  that  are  beyond  the  scope  of  this  study. 

A.  PRIMARY  AND  SECONDARY  BONE 

Primary  bone  is  the  initial  framework  created  for  the  skeletal  system.  Primary 
bone  is  created  from  the  epiphyseal  cartilage  found  near  the  epiphyseal  plate,  also  known 
as  the  growth  plate  [68].  The  mineralization  of  growing  bone  occurs  rapidly,  and  less 
orderly  than  in  secondary  bone.  This  results  in  very  little  correlation  between  collagen 
structure  and  mineral  deposition  [68].  This  disordered  bone  is  very  similar  to  woven  bone 
and  is  generally  stronger  than  secondary  bone  [116].  Secondary  bone  is  formed  when 
primary  bone  is  resorbed  and  new  bone  is  deposited.  Osteoclasts  resorb  bone  while 
osteoblasts  deposit  bone  [71].  Deposition  of  bone  is  directed  by  the  stresses  and  loading 
applied  to  the  bone  [116].  This  occurs  at  a  cellular  level  through  mechanosensation, 
which  allows  the  production  of  new  bone  to  match  the  local  mechanical  environment 
[64],  These  processes  occur  naturally  to  remove  old  bone  and  to  repair  injuries. 

There  are  two  responses  to  bone  injury:  repair  and  regeneration  [117].  Initial 
repair  is  completed  through  quickly  grown  woven  bone.  This  acts  as  a  temporary  solution 
until  regeneration  occurs.  Regeneration  serves  to  replace  the  injured  bone  with  secondary 
bone,  effectively  returning  it  to  its  original  state.  However,  there  are  limits  to  bone 
regeneration.  Wounds  that  are  too  large  may  never  heal  [117].  It  is  in  these  scenarios 
when  medical  solutions  can  aid  in  recovery. 

B.  SYNTHETIC  TISSUE  GROWTH 

The  ability  to  create  a  lab  grown  biomimetic  bone  substitute  has  profound  medical 
implications.  Since  bone  is  a  regenerative  material,  most  injuries  will  heal  themselves. 
However,  the  healing  process  is  very  long  and  requires  setting  the  bones  in  the  proper 
positions.  Current  solutions  to  broken  bones  include  metal  rods,  plates  and  screws  that 
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secure  the  broken  bones  in  the  proper  locations.  The  process  can  be  accelerated  through 
bioresorbable  scaffolds  [118].  These  scaffolds  are  made  of  organic  and  inorganic 
materials  that  can  be  broken  down  by  the  body’s  osteoclasts.  Once  the  scaffolding  is 
resorbed,  new  bone  is  laid  down. 

Research  is  currently  being  done  into  the  feasibility  of  synthetic  bone  grafts,  yet, 
reproducing  the  hierarchical  structure  of  bone  has  proved  difficult.  The  most  difficult 
hierarchy  to  achieve  is  intrafibrillar  mineralization  [68].  The  process  is  directed  by 
proteins  within  the  gap  zone  that  are  not  present  in  laboratory  conditions.  The  lack  of  cell 
regulated  growth  has  limited  laboratory  results. 

There  are  several  factors  influencing  synthetic  bone  composition.  Three  major 
goals  of  scaffold  material  are  biocompatibility,  biodegradability  and  osteoconductivity 
[119].  Biocompatibility  allows  the  scaffold  to  interact  with  the  surrounding  bone, 
biodegradability  allows  it  to  be  broken  down  by  the  osteoclasts  and  osteoconductivity 
allows  bone  structure  to  form  properly.  These  properties  are  highly  depended  on  the 
nano-structure  of  the  synthetic  bone  [119].  In  order  to  mimic  natural  bone  structure  and 
strength,  the  properties  of  bone  at  every  level  should  be  understood. 

C.  AGE  AND  DISEASE 

Both  age  and  disease  can  affect  the  structure  of  bone,  the  effects  of  which  are 
often  similar.  The  most  common  bone  disease  is  osteoporosis.  Osteoporosis  is  a  disease 
that  results  in  reduced  bone  mass  and  density.  This  reduction  of  bone  mass  and  density 
has  a  greater  impact  on  trabecular  bone  than  on  cortical  bone.  The  highly  porous  structure 
of  trabecular  bone  leads  it  to  be  more  metabolically  active.  In  bone  disorders,  this  activity 
leads  to  greater  absorption  rate  than  deposition  rate,  effectively  reducing  bone  mass.  This 
can  result  in  both  trabecular  thinning  and  trabecular  loss  [120],  Trabecular  thinning 
results  in  reduced  radii  of  trabeculae  while  trabecular  loss  is  the  complete  perforation  of 
trabeculae  that  results  in  the  loss  of  a  stmt.  Of  these  two  mechanisms,  trabecular  loss  has 
a  more  profound  impact  on  mechanical  strength. 

Age  and  disease  affect  individuals  differently  and  little  is  known  about  the  cause 
of  age-related  bone  loss  [89].  General  trends  linked  to  gender  are  observed  with  bone 
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loss.  Bone  loss  in  females  is  linked  to  a  decrease  in  estrogen;  the  decrease  of  estrogen 
associated  with  menopause  increases  osteoclast  activity  [89].  This  increases  absorption  in 
relation  to  deposition,  causing  bone  loss.  In  males  the  onset  of  bone  loss  is  linked  to  a 
decrease  of  androgen  [89].  The  onset  of  bone  loss  is  less  rapid  in  males  than  in  females. 
In  addition  to  trends  with  the  onset  of  bone  loss,  there  exist  trends  with  the  physical 
manifestation  of  bone  loss.  Females  experience  greater  trabecular  loss  while  males  are 
prone  to  trabecular  thinning  [89].  Furthermore,  the  onset  of  bone  loss  begins  at  different 
ages  for  males  and  females.  The  process  generally  starts  as  a  rapid  loss  in  the  fifties  for 
women  whereas  men  have  gradual  onset  in  their  sixties  [89]. 

The  model  presented  in  this  study  analyzed  the  effects  of  decreased  bone  density 
on  the  macroscale  properties  of  cancellous  bone.  A  complex  FEM  could  more  accurately 
model  the  properties  of  cancellous  bone.  Not  only  could  the  effects  of  both  trabecular 
thinning  and  trabecular  loss  be  accounted  for,  but  also  the  bi-modulus  property  of 
trabecular  bone. 

D.  IN  VIVO  WATER 

The  presence  of  water  influences  the  structure  and  function  of  living  materials, 
both  in  their  formation  and  their  in  vivo  state.  Mineralization  in  vivo  occurs  in  an 
aqueous  environment,  at  atmospheric  pressures,  and  at  ambient  temperatures  [121]. 
Mineralization  occurs  after  the  fibril  has  formed  its  structure.  The  water  is  tightly  packed 
in  a  dense  liquid  phase  between  molecules  and  bound  within  the  gap  region.  BSP  is 
bound  within  the  dense  liquid  phase  near  the  gap  region  and  nucleates  crystal  growth. 
The  dense  packing  of  water  inhibits  final  crystal  growth  but  promotes  migration  of 
crystal  nucleation.  This  causes  the  HA  crystals  to  propagate  within  the  gap  region  and  to 
form  small,  thin,  plates  aligned  with  the  fibril  [121]. 

The  presence  of  water  affects  the  properties  of  biomaterials  as  well  as  their 
formation.  Variations  in  tissue  hydration  have  a  profound  impact  on  mechanical 
properties.  This  effect  is  due  in  part  to  the  structural  importance  of  water.  Water  helps  to 
bind  the  TC  molecules,  providing  structural  support  through  hydrogen  bonds  [122],  The 
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geometric  relationship  of  TC  and  HA  accounts  for  the  material  properties  of  bone. 
Abnormalities  in  geometry  can  bring  about  changes  in  the  resulting  properties. 

Furthermore,  the  freezing  of  biomaterials  will  have  an  effect  on  their  water  phase. 
This  is  relevant  to  posthumous  testing  of  human  samples,  as  almost  all  samples  are  frozen 
for  storage.  This  topic  deserves  further  research. 

E.  SOFT  TISSUES 

Biomaterials  are  not  limited  to  bone,  as  other  soft  tissues  can  be  readily  identified. 
Among  these,  tendons,  ligaments,  and  muscles  are  the  most  prevalent. 

1.  Tendon  and  Ligaments 

Tendons  and  ligaments  are  a  form  of  connective  tissue.  Tendons  connect  bone  to 
muscle  while  ligaments  connect  bone  to  bone.  Tendons  and  ligaments  are  responsible  for 
the  transmission  of  tensile  loads  required  for  motion  in  the  human  body,  the  restriction  of 
skeletal  movement,  and  the  storage  of  elastic  energy  [123].  Both  tendons  and  ligaments 
are  hierarchical  biomaterials. 

The  largest  structures  of  ligaments  and  tendons  are  fascicles.  Fascicles  are  groups 
of  fibers  surrounded  by  an  inter-fascicular  matrix  [124],  Fibers  are  themselves  a  slightly 
twisting,  lattice  array  of  fibrils  [123].  Fibrils  are  an  ordered  arrangement  of  collagen 
molecules.  Fibrils  exhibit  a  crimp,  which  is  due  to  the  relaxed  state  of  the  collagen 
molecules.  This  crimp  unfolds  when  a  tensile  load  is  applied  to  the  fibril  and  reforms 
when  the  structure  is  relaxed  [124],  The  wavy  crimps  of  the  fibrils  are  aligned,  which 
results  in  a  visible  waveform  within  the  fibers  [124],  Fibrils  are  connected  through 
intermolecular  cross-linking,  as  with  bone  fibrils  [125].  This  cross-linking  produces  a 
periodic  banding.  The  fibrils  of  tendons  are  less  mineralized  than  bone,  resulting  in 
greater  flexibility  and  reduced  rigidity  [5],  Collagen  molecules  of  tendons  and  ligaments 
are  similar  to  TC  molecules  present  in  bone,  as  collagen  I  is  the  predominant  form  of 
collagen.  However,  the  collagen  molecules  of  tendons  and  ligaments  also  contain  small 
amounts  of  collagen  II  and  collagen  III  [124],  The  hierarchy  of  tendons  and  ligaments 
can  be  seen  in  Figure  2 1 . 


82 


Diameter  1 .5  nm 


10  nm 


1 0-350  nm 


1-20  pm 


50-300  pm 


15  mm 


Intrafascicular  /  Nerves  and  / 
tenocytes  /  blood  vessels/ 


Endotenon 


Epitenon 


Figure  21 .  Tendon  and  ligament  hierarchy,  from  [126] 


Tendons  and  ligaments  only  transmit  tensile  loads.  Both  ligaments  and  tendons 
have  a  high  water  content,  which  helps  to  lubricate  individual  fibers  against  longitudinal 
motion  [5].  Additionally,  the  liquid  present  in  tendons  and  ligaments  has  a  higher 
viscosity  than  the  liquid  found  in  bone,  allowing  for  shear  forces  to  be  transmitted  to 
proximal  fascicles  [121].  Furthermore,  both  ligaments  and  tendons  have  limited  blood 
supply  and  nutritional  support,  provided  through  capillary  networks.  This  limited 
capillary  network  reduces  the  void  space  within  tendons  and  ligaments  but  also  slows 
healing. 
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2. 


Muscle 


Muscles  are  the  force  producing  tissues  of  the  human  body.  Three  types  of 
muscles  exist:  cardiac  muscle,  skeletal  muscle  and  smooth  muscle  [71].  Cardiac  muscle 
exists  in  the  heart  and  is  capable  of  synchronized  contraction.  Smooth  muscles  carry  out 
autonomous  motion,  such  as  digestions  and  organ  function.  Skeletal  muscle  is  under 
voluntary  control,  and  is  the  main  motion  inducing  tissue  in  the  human  system.  Muscles 
generate  force  through  contractions  [71],  and  are  built  upon  a  hierarchy  of  structures.  The 
muscle  hierarchies  can  be  seen  in  Figure  22. 
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Figure  22.  Muscle  hierarchy,  from  [127]. 
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The  largest  of  the  structures  is  a  muscle  fascicle.  Fascicles  contain  the  connective 
tissue  that  retains  the  muscle’s  shape,  houses  the  blood  and  nerve  supplies,  and  encircles 
muscle  cells  [71].  Muscle  cells  are  the  next  hierarchy  of  muscles;  they  are  often  called 
muscle  fibers.  Muscle  fibers  are  surrounded  by  a  sarcolemma,  which  contains  muscle 
nuclei  and  mitochondria.  These  help  support  the  metabolic  functions  of  muscle. 
Additionally,  the  sarcolemma  aids  in  signal  transmission  of  the  attached  nerves.  The 
smallest  unit  of  muscle  is  the  myofibril,  which  are  contained  within  the  muscle  fibers. 

Each  myofibril  is  surrounded  by  a  sarcoplasmic  reticulum  that  aids  in  both 
containing  the  myofibril  and  propagating  signals  from  the  nervous  system.  A  myofibril  is 
an  orderly  arrangement  of  thick  and  thin  filaments  [71].  The  filaments  of  muscle  are  the 
structures  that  induce  motion.  Thick  filaments  are  a  twisted  array  of  myosin  molecules. 
These  molecules  have  an  attachment  point  at  their  ends  that  connect  to  the  thin  filaments. 
Thin  filaments  alternate  with  thick  filaments  and  are  made  of  actin  molecules.  The  actin 
molecules  have  a  myosin  binding  site,  which  allows  the  thick  filaments  to  bond  to  the 
thin  filaments.  Each  thick  and  thin  filament  is  oriented  parallel  to  the  longitudinal  axis  of 
the  muscle  and  are  symmetrical  across  their  centers  [71]. 

When  a  signal  from  the  nervous  system  initiates  motion,  the  thick  filaments  bond 
to  the  thin  filaments  and  pull  them  inward.  This  shortens  the  myofibril  length,  inducing  a 
tensile  force.  This  action  is  repeated  for  all  myofibrils,  creating  motion. 


85 


THIS  PAGE  INTENTIONALLY  LEFT  BLANK 


86 


VII.  CONCLUSIONS 


The  structures  of  biomaterials  are  highly  dependent  on  complex  geometries.  This 
study  has  shown  that  material  properties  of  hierarchical  structures  can  be  found  by 
analyzing  each  level  independently.  By  linking  all  hierarchies  and  adjusting  parameters, 
the  influence  of  each  level  can  then  be  analyzed.  Additionally,  changes  to  the  geometry  at 
each  level  can  be  completed  to  test  assumptions  about  the  structure  of  bone. 

The  nanoscale  constituents  of  bone  were  described  and  a  spring  model  was 
utilized  to  calculate  the  longitudinal  stiffness  of  tropocollagen.  This  simple  model 
produced  accurate  results.  Additionally,  a  unit-cell  based  micromechanics  model  was 
used  to  analyze  the  microscale  components  of  bone.  Bone  fibrils  represent  a  particle 
reinforced  matrix  while  bone  fibers  are  a  fiber  reinforced  composite.  The  fibrillar  model 
produced  accurate  results  as  compared  to  accepted  values.  The  fiber  model  could  not  be 
compared  to  experimental  results,  so  the  eight  variations  tested  were  carried  over  to  the 
next  hierarchy. 

The  first  macroscale  structures  of  bone  are  the  lamellar  layers.  These  were 
modeled  as  a  fiber  reinforced  composite  and  compared  to  accepted  values  of  mineral 
fraction  volume.  The  results  disproved  the  twisting  hydroxyapatite  fibrillar  model.  The 
linear  fibril  model  was  utilized  for  both  cortical  and  trabecular  bone.  These  macroscale 
structures  utilized  a  layered  composite  model  to  calculate  their  transverse  isotropic 
material  properties.  The  original  hierarchy  calculated  values  below  what  is  expected. 
Cortical  bone  had  a  wider  gap  between  the  model  results  and  accepted  results  than 
trabecular  bone. 

A  simple  finite  element  model  of  a  tetrakaidecahedron  was  used  to  model 
cancellous  bone.  Utilizing  indices  found  through  experimental  testing,  the  volume 
fractions  of  bone  and  trabecular  thicknesses  were  defined  for  different  anatomical 
locations.  The  model  was  shown  to  accurately  predict  the  properties  of  cancellous  bone. 

This  study  has  revealed  that  the  current  anatomical  knowledge  of  bone 
composition  and  structure  is  lacking  in  detail.  The  control  model  analyzed  the  properties 
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of  bone  based  on  the  current  understanding  of  bone  structure.  But  this  resulted  in  a  low 
mineral  volume  fraction  and  low  stiffnesses.  By  altering  the  parameters  of  each  hierarchy 
an  optimized  material  was  found  that  has  material  properties  close  to  the  accepted  values. 

The  current  model  contains  HA  within  the  gap  regions  of  the  collagen  framework 
and  on  the  exterior  of  fibrils  as  extra  fibrillar  matrix.  However,  this  study  has  shown  that 
more  HA  is  present  in  the  microstractures  of  bone.  Intrafibrillar  mineralization  may  exist 
between  the  TC  molecules  of  fibrils.  This  would  stiffen  the  lowest  microscale  level, 
which  is  then  compounded  throughout  the  macroscale  hierarchy.  Alternatively,  the 
extrafibrillar  matrix  may  take  up  a  larger  volume  percent  than  believed.  This  would 
increase  the  mineral  volume  fraction  of  bone,  as  shown  through  Model  4  of  the  analysis. 

Additionally,  it  has  been  assumed  that  the  structure  of  cortical  bone  and  trabecular 
bone  does  not  vary  at  the  lower  hierarchies  [9].  Current  beliefs  attribute  the  difference  in 
material  properties  to  the  varying  orientation  of  lamellar  layers.  But,  this  model  has 
disproved  that  assumption.  The  accepted  results  of  cortical  bone  show  an  average  of  17 
GPa  and  21  GPa  for  the  longitudinal  and  transverse  stiffness,  respectively.  The  accepted 
results  of  trabecular  bone  exhibit  approximately  14  GPa  for  transverse  and  longitudinal 
stiffness.  The  results  of  the  optimized  hierarchy  show  that  the  upper  limit  of  cortical  bone 
stiffness  is  13  GPa  and  15  GPa  for  the  longitudinal  and  transverse  directions, 
respectively.  However,  trabecular  bone  exhibits  17  GPa  and  24  GPa  for  longitudinal  and 
transverse  stiffness,  respectively.  The  only  difference  between  trabecular  and  cortical 
bone,  in  this  study,  is  the  fiber  distribution  and  the  volume  fraction  assigned  to  void 
space.  Yet,  the  values  for  trabecular  bone  are  much  greater  than  cortical  bone.  This 
suggests  a  structural  difference  at  the  lower  hierarchies. 

This  model  can  be  used  to  validate  future  discoveries  about  the  structure  of  bone. 
As  technology  advances,  imaging  capabilities  will  allow  the  nanostructures  of  bone  to  be 
explored  in  detail.  The  discoveries  can  be  checked  against  this  model  to  assess  their 
impact  on  macroscale  properties.  Additionally,  the  properties  of  synthetic  bone  materials 
can  be  checked  against  the  hierarchical  structure  of  bone.  This  would  allow  for  more 
anatomically  beneficial  bone  grafts. 
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APPENDIX.  TETRAKAIDECAHEDRON  FINITE  ELEMENT 
METHOD  STIFFNESS  AND  TRANSFORMATION  MATRICES 


A.  ELEMENT  MATRICES  FOR  THE  THREE  DIMENSIONAL  GLOBAL 
MATRIX,  FROM  [101]. 
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[*22  ] 
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The  variables  used  in  stiffness  matrices  are:  E  which  represents  Young’s  modulus 
along  the  long  axis  of  the  element,  G  which  represents  shear  modulus  of  the  element,  A 
which  represents  cross-sectional  area  of  the  element,  I  which  represents  moment  of 
inertia  in  the  axis  noted  in  subscript,  and  /  which  represents  element  length. 
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B. 


TRANSFORMATION  MATRIX,  FROM  [101]. 
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The  vectors  utilized  in  the  transformation  matrix  are  based  on  the  direction 
cosines  of  the  beam  elements.  Vj  is  directed  along  the  long  axis  of  the  element.  v,is 

directed  perpendicular  to  the  long  axis  of  the  element.  v3  is  perpendicular  to  both  Vj  and 
v, .  The  three  vectors  are  bound  in  three  dimensional  space  by  a  v0  vector,  which  points 
toward  the  origin  of  the  global  axes  from  the  starting  point  of  the  vector.  The  three 
direction  cosine  vectors  uphold  the  following  equations: 

v2  =V1XV0 
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